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1  INTRODUCTION 


In  numerical  antenna  modeling,  the  charge  and  current  distribution  can  have  con¬ 
vergence  problems.  The  charge  and  current  distribution  may  not  converge  as  the  number  of 
basis  functions  in  the  model  are  increased.  Dudley  (1985)  discusses  mathematical  and 
numerical  problems  in  the  method-of-moments  solution  of  antenna  problems.  The  electro¬ 
magnetic  method-of-moments  approach  has  mathematical  limitations  that  contribute  to 
the  convergence  problems.  These  mathematical  problems  are  not  present  in  the  electro¬ 
static  method-of-moments  approach.  Insight  into  the  limitations  of  the  method-of- 
moments  approach  can  be  gained  by  considering  electrically  small  antennas.  The  algo¬ 
rithm  developed  in  this  report  improves  the  convergence  of  the  method-of-moments 
approach.  This  algorithm  allows  the  error  in  the  boundary  condition  to  be  estimated  on 
each  segment. 

The  numerical  limitations  of  the  method  of  moments  can  be  understood  by  consider¬ 
ing  electrically  small  antennas.  Electrically  small  antennas  can  be  analyzed  with  the  quasi¬ 
static  approach  or  with  the  full-wave  approach.  In  principle,  the  full-wave  approach  should 
have  the  same  accuracy  and  stability  properties  as  the  quasi-static  approach.  The  only 
difference  in  the  two  methods  is  the  numerical  formulation  of  the  problem. 

In  the  quasi-static  approach,  the  interaction  between  the  charges  is  divided  into  two 
parts:  the  static  and  frequency-dependent  interaction  terms.  The  static  interaction  terms 
are  much  larger  than  the  frequency-dependent  interaction  terms.  The  static  interaction 
term  is  used  to  compute  the  static-charge  distribution.  The  current  and  frequency- 
dependent  part  of  the  interaction  give  a  small  correction  to  the  static-charge  distribution. 
The  quasi-static  approach  is  numerically  very  stable  and  accurate;  however,  it  is  limited  to 
electrically  small  antennas. 

In  the  full-wave  method-of-moments  approach,  the  static  and  frequency-dependent 
interaction  terms  are  combined  in  one  matrix,  the  impedance  matrix.  For  electrically 
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small  antennas,  the  frequency-dependent  terms  cancel  out  in  the  solution  of  the  impedance 
matrix.  The  cancellation  of  the  frequency-dependent  terms  leads  to  round-off  errors  in  the 
solution  of  the  matrix.  The  full-wave  method-of-moments  approach  is  not  as  stable  or 
accurate  as  the  quasi-static  approach. 

The  basis  functions  and  boundary  condition  developed  in  this  report  simplify  the 
structure  of  the  impedance  matrix.  For  electrically  small  antennas,  the  impedance  matrix 
has  large  static  interaction  terms  with  small  frequency-correction  terms.  The  quasi-static 
properties  of  electrically  small  antennas  are  explicit  in  the  structure  of  the  impedance 
matrix.  This  improves  the  stability  and  accuracy  of  the  full-wave  approach  without  limit¬ 
ing  the  scope  of  the  applications.  This  approach  provides  a  better  understanding  of  the 
error  in  the  method  of  moments. 

The  idea  for  the  numerical  algorithm  described  in  this  report  developed  out  of  an 
effort  to  improve  the  convergence  of  MININEC  3.  The  MININEC  3  code,  developed  by 
Logan  and  Rockway  (1986),  uses  square  pulses  to  model  antenna-current  distributions. 
The  square-pulse  charge  distributions  are  implicit  to  the  MININEC  3  code.  This  simple 
model  gives  surprisingly  good  results.  The  problem  with  a  square-pulse  charge  distribution 
is  that  it  does  not  eliminate  the  parallel  component  of  the  electric  field  on  individual  seg¬ 
ments  of  the  wire.  In  nature,  the  charges  would  redistribute  themselves  to  eliminate  the 
parallel  component  of  the  electric  fields. 

The  charge  and  current  basis  functions  discussed  in  this  report  emulate  nature.  The 
net  charge  on  a  segment  is  modeled  by  a  square  pulse.  The  higher  order  basis  functions  are 
used  to  redistribute  the  charge  on  each  segment  to  eliminate  the  electric  field  along  the 
segments.  Legendre  polynomials'*'  (Richmond,  1965)  are  used  for  the  charge  basis  func¬ 
tions.  The  above  observations  were  the  motivation  for  the  numerical  model  developed  in 


1.  J.  H.  Richmond  expanded  a  current  distribution,  computed  with  sine  and  cosine  basis 
functions,  with  Chebyshev  and  Legendre  polynomial  series,  which  give  faster  convergence. 
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this  report.  The  general  antenna  problem  is  too  complex  for  an  initial  numerical  electro¬ 
magnetic  code.  A  numerical  code  for  the  symmetric  dipole  is  used  to  test  this  approach  to 
numerical  electromagnetic  modeling. 

In  section  2,  the  stability  and  convergence  of  this  algorithm  are  illustrated  with  a 
computer  model  of  a  dipole.  The  admittance  at  resonance  and  antiresonance  of  the  dipole 
is  computed  using  one  to  seven  basis  functions  per  segment.  The  admittance  for  the  one, 
two,  and  seven  basis  function  models  are  plotted  against  the  matrix  size.  The  admittance 
curves  agree  with  King’s  (1956)  calculations  for  the  resonance  and  antiresonance  of  the 
dipole.  Their  convergence  and  accuracy  are  a  consequence  of  the  new  boundary  condition. 

For  electrically  small  antennas,  the  charge  distribution  can  be  compared  to  the 
static  charge  distribution  (Watt,  1967).  A  24x24  and  a  48x48  impedance  matrix  were  used 
to  show  the  convergence  of  the  charge  distribution  model.  Different  combinations  of  seg¬ 
ments  and  basis  functions  were  used  to  model  the  charge  distribution.  A  simple  model 
with  two  basis  functions  per  segment  gives  an  accurate  model  of  the  static  charge  distribu¬ 
tion  on  the  interior  segments  in  the  model.  A  model  with  eight  basis  functions  per  segment 
gives  the  best  results  for  the  charge  distribution  on  the  entire  antenna.  The  accuracy  of 
the  model  at  the  end  of  the  wire  is  improved  by  use  of  a  84x84  impedance  matrix  with 
12  segments  and  seven  basis  functions  per  segment.  These  results  show  that  the  charge 
distribution  converges  as  the  number  of  basis  functions  is  increased.  The  numerical  results 
indicate  a  very  stable  model  of  the  thin-wire  dipole. 

The  development  of  the  algorithm  includes  a  discussion  of  the  different  types  of 
errors  in  the  algorithm.  The  errors  fall  into  two  classes:  modeling  and  computational 
errors.  The  basis  functions  used  to  model  the  charge  distribution  are  only  an  approxima¬ 
tion  of  the  correct  charge  distribution.  The  second  modeling  error  is  caused  by  the  approx¬ 
imation  of  the  ideal  boundary  condition.  Error  in  the  boundary  condition  induces  errors  in 
the  charge  and  current  distributions.  The  interaction  of  the  two  modeling  errors  is  deter¬ 
mined  by  the  numerical  algorithm  used  in  computing  and  imposing  the  boundary  condi- 
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tion.  The  numerical  algorithm  includes  several  sources  of  computational  errors:  evaluation 
of  the  fields,  imposing  the  boundary  condition,  and  solution  of  the  impedance  matrix.  The 
following  discussion  is  a  qualitative  discussion  of  the  interaction  of  the  different  errors  in 
our  numerical  formulation  of  the  dipole  antenna  model. 

The  notation  is  developed  in  section  3.  This  notation  allows  a  detailed  discussion  of 
the  charge  and  current  distribution  on  each  segment  in  the  antenna.  The  scalar  and  vector 
potentials  on  each  segment  are  a  function  of  the  charge  and  current  basis  function  on  the 
other  segments.  The  scalar  and  vector  potentials  are  used  to  formulate  the  boundary  con¬ 
dition  on  each  segment. 

The  Legendre-polynomial  charge  basis  functions  are  discussed  in  section  4.  The 
current  basis  functions  are  computed  from  the  continuity  equation  and  the  charge  basis 
functions.  The  charge  distribution  is  modeled  by  N  basis  functions,  Pg(t)  to  P^_^(t),  for 
each  segment  of  the  antenna.  The  charge  distribution  modeling  error  is  PN^)’  *^+1^’ 
plus  additional  higher  order  basis  functions.  The  modeling  error  in  the  charge  distribution 
causes  errors  in  the  scalar  and  vector  potentials. 

Section  5  shows  that  the  Legendre  polynomials  are  a  multipole  expansion  of  the 
charge  distribution  on  the  segment.  The  potential  integral  of  each  basis  function  contrib¬ 
utes  to  the  far,  intermediate,  and  near  fields.  The  potential  in  each  region  is  discussed  in 
detail.  The  static  part  of  the  potential  falls  off  as  r  for  the  nfch  Legendre  polyno¬ 

mial.  The  modeling  error  in  the  charge  distribution  impacts  the  potentials  and  boundary 
condition  only  on  nearby  segments.  This  algorithm  is  compared  to  the  matrix  transforma¬ 
tion  technique  developed  by  Canning  (1992).  The  numerical  calculation  of  the  scalar  and 
vector  potentials  is  discussed  in  detail  in  section  6. 

The  scalar  and  vector  potentials  are  computed  from  the  Legendre  polynomial  poten¬ 
tial  integral,  which  is  the  integral  of  the  surface  charge  distribution  and  the  free-space 
Green’s  function.  The  integral  is  evaluated  by  integrating  along  the  length  of  the  segment 
first,  and  then  integrating  around  the  circumference.  The  integration  around  the  circum- 
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ference  is  evaluated  by  Gauss-Chebyshev  integration,  which  has  a  very  simple  physical 
interpretation.  It  is  equivalent  to  replacing  the  surface  charge  with  equally  spaced  lines  of 
charge  parallel  to  the  axis  of  the  cylinder.  This  simple  physical  interpretation  of  the  poten¬ 
tial  integral  is  used  to  deduce  the  error  term  for  the  numerical  integration. 

A  sequence  of  Gauss-Chebyshev  integration  formulas  can  be  used  to  estimate  the 
error  in  the  angular  integration.  Some  examples  are  used  to  illustrate  the  convergence  of 
the  numerical  evaluation.  In  addition,  the  static  potential  near  the  end  of  the  segment  is 
plotted  to  show  the  sharp  drop  in  the  potential  at  the  end  of  the  segment.  The  potential 
integral  can  be  computed  with  arbitrary  accuracy.  The  calculation  of  the  potential  integral 
does  not  introduce  errors  into  the  boundary  condition. 

In  principle,  the  electric  field  should  be  normal  to  the  surface  of  a  perfect  conductor 
at  all  points  on  the  surface.  In  the  method  of  moments,  the  ideal  boundary  is  replaced  by 
an  approximation  of  the  boundary  condition.  The  modeling  errors  introduced  by  the 
boundary  condition  are  discussed  in  section  7.  In  addition,  we  discuss  numerical  errors 
introduced  by  the  numerical  evaluation  of  the  boundary  condition.  The  numerical  and 
modeling  errors  primarily  create  high  multipole  errors  in  the  charge  distribution.  These 
errors  impact  only  the  nearby  segments.  The  numerical  method  used  to  calculate  the 
boundary  condition  can  be  used  to  estimate  the  error  in  the  boundary  condition. 

In  section  8,  we  derive  the  impedance  matrix  created  by  these  basis  functions  and 
the  boundary  conditions.  The  impedance  matrix  neglects  some  high  multipole  moment 
terms  introduced  by  the  current  basis  functions.  The  errors  caused  by  this  approximation 
are  local  to  each  segment.  The  local  error  on  each  segment  can  be  computed  by  extending 
the  impedance  matrix  to  include  extra  boundary  condition  expansion  coefficients.  These 
extra  boundary  condition  expansion  coefficients  represent  the  error  in  the  boundary  condi¬ 
tion  on  each  segment. 

The  stability  of  the  impedance  matrix  is  evaluated  by  computing  the  condition 
number.  The  condition  number  depends  on  the  number  of  segments  used  in  the  model;  it 
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does  not  depend  on  the  number  of  basis  functions  used  on  each  segment.  This  formulation 
of  the  impedance  matrix  is  compared  to  the  polynomial  numerical  code  developed  by 
Djordjevic  et  al.  (1990). 

Section  9  gives  the  summary  and  the  conclusions.  A  list  of  variables  and  their  defi¬ 
nitions  is  given  in  appendix  B. 
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2  NUMERICAL  RESULTS 


The  stability  and  convergence  of  our  numerical  algorithm  is  illustrated  with  the 
numerical  calculations  for  a  dipole.  The  admittance  of  the  dipole  calculated  by  using  our 
algorithm  is  compared  to  the  analytic  approximation  of  King  (1956).  A  4-  and  8-segment 
dipole  models  give  almost  identical  admittances  over  a  wide  range  of  frequencies.  Detailed 
convergence  tests  are  also  run  for  seven  different  models  at  the  resonance  and  antiresonance 
of  the  dipole.  These  results  are  compared  to  King’s  (1956)  results.  The  numerically  calcu¬ 
lated  charge  distribution  is  compared  to  the  static-charge  distribution  for  electrically  small 
dipoles. 

The  dipole  is  described  by  three  parameters:  the  total  length,  2h  =  lm;  the  wire 
radius,  a  =  4.5401E— 5;  and  the  gap  at  the  feed  point,  d  =  0.  In  King’s  (1956)  notation,  Q 
is  20,  where  0  =  2  ln(2h/a).  Each  arm  of  the  dipole  is  divided  into  R  segments.  The  pro¬ 
gram  uses  the  symmetry  of  the  dipole  to  eliminate  the  unknowns  on  one  arm  of  the  dipole. 
The  algorithm  can  use  one  to  seven  basis  functions  per  segment. 

Figure  2—1  plots  the  admittance  for  an  8-segment  model  with  three  basis  functions 
per  segment  and  King’s  (1956)  admittance;  the  matrix  size  is  24x24.  The  results  are  in 
very  close  agreement.  The  convergence  of  the  algorithm  is  illustrated  by  comparing  the 
admittance  computed  with  4-  and  8-segment  models  with  three  basis  functions  per  seg¬ 
ment;  the  matrix  sizes  are  12x12  and  24x24,  respectively.  In  this  model,  a  gap  of  d  =  0.01 
meters  was  used  to  eliminate  the  logarithmic  singularity  in  the  charge  distribution.  A 
constant  current  is  used  on  the  source  segment  in  the  gap.  The  admittance  for  4-  and  8- 
segment  models  is  plotted  in  figure  2-2  as  a  function  of  dipole  length  measured  in  wave¬ 
lengths.  The  4-segment  model  uses  a  total  of  8  segments  to  model  both  arms  of  the  dipole. 
The  numerical  calculation  uses  only  2-segments  per  wavelength  at  a  dipole  length  of  four 
wavelengths.  The  error  in  the  2-segment s-per-wavelength  model  is  very  small.  The  small 
difference  between  the  two  admittance  curves  is  plotted  in  figure  2—3.  The  error  increases 
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Figure  2-3.  Difference  in  admittance  for  the  4-  and  8-segment  models  with  three  Legendre  polynomials 
per  segment. 


as  fewer  segments  are  used  per  wavelength. 

The  convergence  of  several  dipole  models  is  evaluated  at  the  resonance  frequency  of 
the  antenna.  The  resonance  of  a  1-meter  dipole  is  at  146.0  MHz.  The  admittance  is  com¬ 
puted  for  one  to  seven  basis  functions  per  segment.  The  conductance  and  susceptance  are 
plotted  in  figures  2—4  and  2—5  as  functions  of  impedance  matrix  size.  The  size  of  the  impe¬ 
dance  matrix  is  determined  by  the  number  of  unknowns  used  to  model  the  charge  distribu¬ 
tion.  The  plot  shows  the  results  for  three  different  charge  distribution  models:  one,  two, 
and  seven  basis  functions  per  segment.  The  charge  distribution  models  with  three  to  six 
basis  functions  per  segment  are  bound  by  the  two  and  the  seven  basis  function  curves.  The 
two  and  the  seven  basis  function  conductance  curves  are  within  1%  of  King’s  (1956)  com¬ 
puted  result,  the  dotted  line.  The  error  in  the  one  basis  function  conductance  curve  is  less 
than  1%  for  nine  or  more  segments.  The  percentage  error  for  the  susceptance  cannot  be 
computed.  However,  the  error  in  the  susceptance  can  be  viewed  as  a  phase  error  in  the 
admittance.  The  phase  error  is  less  than  1.2°  for  the  two  and  seven  basis  function  suscep¬ 
tance  curves.  The  phase  error  in  the  one  basis  function  susceptance  curve  is  less  than  1° 
for  the  nine  or  more  segments  in  the  model.  The  three  curves  converge  to  the  same  con¬ 
ductance  values,  the  typical  difference  in  the  conductance  curves  is  7/100%  of  King’s 
(1956)  value. 

The  convergence  of  several  dipole  models  is  also  evaluated  at  the  antiresonance 
frequency  of  the  dipole.  The  antiresonance  of  a  1-meter  dipole  is  at  281.51  MHz.  The 
conductance  and  susceptance  are  plotted  in  figures  2—6  and  2—7  as  functions  of  the  impe¬ 
dance  matrix  size.  The  plot  shows  the  results  for  three  different  charge  distribution 
models:  one,  two,  and  seven  basis  functions  per  segment.  The  charge  distribution  models 
with  three  to  six  basis  functions  per  segment  are  bound  by  the  two  and  the  seven  basis 
function  curves.  The  seven  basis  function  conductance  curves  are  within  3/10%  of  King’s 
(1956)  computed  result,  the  dotted  line.  The  two  basis  function  conductance  curves  are 
within  3/10%  of  King’s  (1956)  result  when  two  or  more  segments  are  used  in  the  model. 
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Figure  2-4.  The  three  curves  show  the  convergence  of  the  conductance 
at  resonance  for  one,  two,  and  seven  basis  functions  per  segment. 
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Figure  2-7.  The  three  curves  show  the  convergence  of  the  susceptance  at  antiresonance 
for  one,  two,  and  seven  basis  functions  per  segment. 


The  error  in  the  one  basis  function  conductance  curve  is  less  than  1%  for  six  or  more  seg¬ 
ments.  The  phase  error  is  less  than  2.6°  for  the  one,  two,  and  seven  basis  function  suscep- 
tance  curves  when  14,  3,  and  1  or  more  segments  are  used  in  the  respective  models.  The 
three  curves  converge  to  the  same  conductance  values;  the  typical  difference  in  the  conduc¬ 
tance  curves  is  3/100%  of  King’s  (1956)  value. 

The  above  plots  indicate  that  the  current  converges  at  the  feed  point.  The  current 
at  the  feed  point  determines  the  net  charge  ( ju )  in  each  arm  of  the  antenna;  however,  the 
charge  distribution  could  be  incorrect.  The  convergence  of  the  charge  distribution  is  a  very 
sensitive  test  of  the  stability  of  our  model.  In  our  model,  the  charge  distribution  at  the 
junction  of  two  segments  is  not  required  to  be  continuous.  A  discontinuity  in  the  charge 
distribution  would  create  an  error  in  the  boundary  condition.  The  magnitude  of  the  discon¬ 
tinuity  in  the  charge  distribution  is  a  good  test  of  the  accuracy  of  the  solution.  Several 
different  polynomial  models  were  tested  for  convergence  and  accuracy. 

The  charge  distribution  for  an  electrically  small  dipole  was  given  in  Watt  (1967).  In 
the  following  plots,  the  feed  point  is  at  0  and  the  end  of  the  wire  is  at  0.5.  Figure  2—8 
shows  the  results  for  a  two-basis  function  model  of  the  charge  distribution  with  12-  and  24- 
segments,  for  a  total  of  24  and  48  unknowns,  respectively.  The  charge  distribution  in  the 
interior  segments  is  continuous.  On  the  feed  point  and  end  segment,  the  error  in  the  linear 
charge  distribution  model  can  be  approximated  with  a  quadratic  charge  distribution.  This 
quadratic  charge  distribution  gives  a  quadrupole  error  in  the  potential  surrounding  the 
segment.  The  charge  distribution  error  on  the  feed  point  and  end  segment  does  not  affect 
charge  distribution  on  the  adjacent  segment.  The  agreement  between  the  12-  and  24- 
segment  models  shows  that  the  algorithm  converges  as  the  number  of  segments  is 
increased.  The  accuracy  of  the  feed  point  and  end  segment  is  improved  by  using  high-order 
polynomial  models.  Figure  2—9  shows  results  for  a  3-  and  6-segment  model  with  eight  basis 
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functions  per  segment  ,  for  a  total  of  24  and  48  unknowns,  respectively.  The  error  at  the 
feed  point  and  end  segment  is  much  smaller  in  this  model.  At  the  junction  of  the  seg¬ 
ments,  the  charge  distribution  has  a  small  jump  in  value.  The  small  difference  between  the 
3-  and  6-segment  models  illustrates  how  the  solution  converges  as  the  number  of  segments 
is  increased. 

The  accuracy  of  the  model  at  the  end  of  the  wire  can  be  improved  by  using  more 
segments.  Figure  2—10  shows  the  charge  distribution  at  the  end  of  the  wire;  the  model  uses 
12  segments  with  seven  basis  functions  per  segment  or  a  84x84  matrix.  This  model  is  very 
accurate  near  the  end  of  the  antenna.  The  accuracy  of  the  model  at  the  end  of  the  wire  is 
attributed  to  the  boundary  condition,  which  is  applied  to  the  entire  segment. 

The  plots  of  the  charge  distribution  indicated  a  very  stable  algorithm.  The  errors 
on  the  feed  point  and  end  segment  do  not  affect  neighboring  segments.  This  stability  can 
be  attributed  to  the  boundary  condition  used  in  this  algorithm. 


2.  The  8-segment  model  neglected  the  current  term  for  the  eight  basis  function.  For 
electrically  small  antennas  the  error  is  insignificant. 
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seven  basis  functions  per  segment. 


3  ELECTROMAGNETIC  EQUATIONS  AND  NOTATION 


The  algorithm  solves  for  the  current  and  charge  distributions  caused  by  an  input 
voltage  at  the  feed  point  of  the  antenna.  The  algorithm  is  implemented  for  a  dipole 
antenna  with  a  length,  2h,  and  a  wire  radius,  a.  The  dipole  is  fed  at  the  center  with  a 
single-frequency  voltage  source.  The  current  flow  at  the  feed  point  creates  a  buildup  of 
charge  on  each  arm  of  the  antenna.  The  charge  on  the  antenna  moves  in  response  to  the 
parallel  component  of  the  electric  field  on  the  surface  of  the  wire.  On  a  perfect  conductor, 
the  charge  moves  to  eliminate  the  parallel  component  of  the  electric  field.  The  algorithm 
developed  in  this  report  could  be  used  to  model  symmetric  thin-wire  antennas  without 
closed  loops. 

The  continuous  charge  distribution  on  the  dipole  is  modeled  in  the  computer  with  a 
finite  number  of  basis  functions.  Each  arm  of  the  dipole  is  divided  into  R  segments  with  a 
length,  L  =  h /R  (figure  3—1).  Segment  1  is  at  the  feed  point,  and  segment  R  is  at  the  end 
of  the  wire.  The  position  on  each  segment,  5,  is  parameterized  by  t  ,  where  — 1  <  t  <  1,  t 
increases  in  the  direction  of  increasing  s.  The  charge  and  vector  current  distributions  on 
segment  s  are  modeled  with  N  basis  functions,  Q  (t  )  and  J  (t  ),  respectively.  The 
italic  subscript  s  indicates  the  segment  number,  and  the  latin  subscript  n  identifies  the 
basis  function  number.  Bold  italic  variables  indicate  a  vector  quantity.  The  unit  vector, 
u  ,  for  each  segment,  lies  on  the  axis  of  the  segment  and  points  in  the  direction  of  increas- 
ing  t^  .  The  current  basis  function,  Jg  g(t ^),  overlaps  with  the  adjacent  segment,  b  =  s—  1. 
In  terms  of  the  basis  functions,  the  charge  and  current  distributions  are  approximated  as 
follows: 


n  =  0 


(3-1) 
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(3-2) 


N— 1 

J( t,)~  £  6  J  (t,)  for  b  =  s  —  1  and  s. 
sK  b  '  sn  snl  o' 

n  =0 


The  current  is  zero  at  the  free  end  of  each  wire  in  the  antenna.  Equation  (3—2)  must  van¬ 
ish  at  the  end  of  the  wire.  The  current  coefficients,  S  ,  and  charge  coefficients,  r?  ,  are 
related  by 


V  •  J 

j  oj 


(3-3) 


The  continuity  equation  can  be  used  to  reduce  2N  unknowns  to  N  unknowns  for  each  seg¬ 
ment  of  the  thin  wire.  The  symmetry  of  the  charge  and  current  distributions  is  used  to 
eliminate  the  unknowns  on  the  left  arm  of  the  dipole. 

The  scalar  potential,  <j>5  n(r),  and  vector  potential,  Ag  n(r),  for  each  basis  function  is 
computed  with  the  free-space  Green’s  function: 


^5  n^r)  Tt? 


47rer 
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-jk\r-r' 


dV', 


r—r' 


A-{r)  =  ±>e 
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J  (r'V 
s  nv  ’  i 


-jk\r-T' 


dV', 


V' 


r—r' 


(3-4) 


(3-5) 


where  ris  the  position  of  the  field  point,  r'  is  the  position  of  the  source,  |  r—r'  |  is  the  dis¬ 
tance  from  source  to  the  field  point,  and  k  is  the  wave  number  of  the  exciting  signal.  The 
free-space  Green’s  function  is  a  point  source  solution  of  Maxwell’s  equations.  The  Green’s 
function  is  integrated  over  the  region  containing  the  source  of  the  field  r',  p(r'),  and  J{r'). 
The  units  are  chosen  as  rationalized  (MKS),  in  free  space,  l/(47re0)  is  10  c  . 
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The  scalar  and  vector  potential  solutions  to  Maxwell’s  equations  are 


and 


V2A  +  A  = 

o  2 


(3-6) 


(3-7) 


The  Lorentz  gauge  condition  is 


c2V-A  —  jhxj)  =  0. 


(3-8) 


The  Lorentz  gauge  condition  has  been  used  to  simplify  the  current  term  (Feynman, 
Leighton  and  Sands,  1964)  on  the  right  side  of  equation  (3—7).  Equations  (3—6)  and  (3—7) 
relate  the  potentials  to  the  general  charge  and  current  distributions,  p  and  J,  where  p  and  J 
obey  the  continuity  equation  (3—3). 

The  charge  and  current  distributions,  p  and  J,  are  determined  indirectly  with  an 
electric-field  boundary  condition.  The  electric  field  for  perfect  conductors  must  be  normal 
to  the  conducting  surface.  On  each  arm  of  the  antenna,  the  parallel  component  of  the  elec¬ 
tric  field  is  zero  on  the  surface.  At  the  feed  point,  the  input  voltage  is  represented  as  a 
parallel  component  of  the  electric  field. 

The  error  in  the  boundary  condition  is  the  line  integral  of  the  electric  field  from  the 
center  of  the  feed  point,  z  =  0,  to  the  point  Xj  on  segment  /: 


EMF(ty)  = 


•*(ty)  z(*f) 

E-dz  =  <J>(0)  —  <j>(ty)  —  ju  A(z)-dz. 


J0 


J0 


(3-9) 
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The  subscript  /  denotes  the  field  point  on  the  field  segment.  The  potentials,  <j)(ty)  and 
A(z),  are  calculated  with  equations  (3^4)  and  (3—5).  The  approximations  used  in  the 
charge  and  current  basis  functions  introduce  errors  in  the  EMF(ty): 

R  N-l 
EMF(t  r)  =  S  £ 

J  s  =  1  n  =o 

(3—10) 
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s  n 


ju 


Z{tf) 

6  A  (z)'dz 
s  n  sn'  ' 


Equation  (3—10)  relates  the  current  and  charge  distributions  on  the  source  segment 
s  to  the  electromotive  force  (EMF)  on  field  segment  /.  Ideally,  this  integral  is  a  constant 
equal  to  half  of  the  feed-point  voltage  for  a  dipole.  In  practice,  the  continuous  boundary 
condition  is  approximated  by  a  discrete  model. 

The  boundary  condition  on  each  segment  is  imposed  by  expanding  the  EMF(t  A 
equation  (3—10),  on  segment  /as  a  sum  of  Legendre  polynomials 

N-l 

EMF(t,)  =  S  a,  P  (tA  (3-11) 

J  m=0  J  J 


The  N  boundary  conditions  for  each  segment  are  imposed  by  setting  the  expansion  con¬ 
stants  equal  to 
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and 


=0  for  0  <  m  <  N. 

/m 


(3-12) 
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In  equation  (3—13),  the  integral  on  the  right  is  set  equal  to  the  constant,  aym,  defined  in 
equation  (3—12): 


afm  =  <m+5) 


,1 


Pm(t/)EMF(tf)dt 


(3-13) 


Equation  (3—13)  represents  N  boundary  conditions  for  segment  /.  The  total  number  of 
boundary  conditions  for  the  antenna  is  R-  N. 

Equation  (3—10)  is  substituted  into  equation  (3—13)  to  give  an  R-N*R-N  interac¬ 
tion  matrix: 


„1 

•2(t f  ) 

“/m  = 

R  N-l 

r  i 

s  =  1  n  =  0 

vw-wvi 

8  A  (z)-dz 
s  n  snv  ; 

-1 

* 

0 

(3-14) 


The  continuity  equation  (3—3)  is  used  to  eliminate  either  Sg  n  or  Vs  n  from  the  above  equa¬ 
tion.  Each  boundary  condition  is  a  function  of  the  N'basis  function  coefficients,  6  n  or 

^5  n‘ 

Equation  (3-14)  is  a  matrix  equation: 


afm  Zfs  mn^-s  n’  ^  ^ 


where  y  =  8  or  ri  The  elements  in  the  matrix,  Z  f  _  ,  will  have  units  of  impe- 
As  n  s  n  's  n  P  mn 

dance  when  y  =8  and  capacity  when  y  „  =  rj  .  The  total  number  of  unknowns  is 

As  n  5  n  J  As  n  's  n 

jR-N.  A  combination  of  8  and  7?  is  used  to  model  the  current  and  charge  distribution. 

s  n  's  n  ° 

For  simplicity,  the  complete  matrix  will  be  called  the  impedance  matrix. 
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4  BASIS  FUNCTIONS 


Harrington  (1968)  was  the  first  to  formulate  the  electromagnetic  problem  in  the 
method  of  moments.  In  his  approach,  each  wire  in  the  antenna  is  divided  into  several  seg¬ 
ments  and  the  charge  and  current  distributions  on  each  segment  of  the  antenna  are  mod¬ 
eled  with  square  pulses.  As  a  result,  both  the  current  and  charge  distributions  are  discon¬ 
tinuous.  The  current  model  can  be  made  continuous  by  introducing  triangular  current 
pulses;  however,  the  charge  distribution  remains  discontinuous.  A  smooth  charge  distribu¬ 
tion  can  be  modeled  by  introducing  extra  current  and  charge  basis  functions.  These  basis 
functions  provide  a  more  realistic  model  of  the  charge  distribution  on  a  wire. 

The  first  current  basis  function  is  a  triangular  basis  function  that  gives  a  square- 
pulse  charge  distribution.  The  triangular  basis  function  shown  in  figure  4— l(row  b)  has 
been  used  by  many  authors  (Michalski  and  Butler,  1983;  Hwu  and  Wilton,  1988).  The 
triangular  basis  function,  Jg  Q(ty),  overlaps  the  two  segments  s-1  and  s: 

Js  o  (*8-1 ) =  ~  +  Vi)  Vu  (4-la) 

Js  0^5  )  =  ~  “  *8  (4-lb) 


J5  o(t/)  =  0  for  ft  S~1  or  s-  (4— lc) 

The  continuity  equation  (3—3)  is  used  to  compute  the  charge  distribution,  Q5  j(t^  ) 


Qs_i  j  {%s-i )  =  !/L>  (4— 2a) 

Q5j(tS)  =  -1/L’  (4"2b) 
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Q$J{tf)  =  0  for  /#  5  or  5-  1. 


(4— 2c) 


The  charge  distribution  described  by  equation  (4—2)  is  a  portion  of  the  charge  distri¬ 
bution  on  the  segment  s.  The  complete  model  of  the  charge  distribution  on  the  segment  s 
uses  two  triangular  pulses,  Jg  g(t^  )  and  ^  g(t5 ).  The  total  charge  on  the  segment, 
Qs(ts),  is 


5J±1_0  {s0 

L 


(4-3) 


which  is  equivalent  to 


“  ^sO^sO^s) 


(4-4a) 


where 


Vs  0  = 


6s+l  0  Ss  0 


and  Q5  0(t5  )  =  PQ(ts  )  for  -1  <  t,  <  1- 


(4-4b) 


Pn(t  )  is  the  first  Legendre  polynomial  with  the  standard  parameterization,  P  (1)  =  1. 
The  constant,  rjs  q,  is  the  derivative  of  the  triangular  current  basis  functions. 

The  charge  and  current  distribution  created  by  the  first  current  basis  function, 
Js  g(t5 ),  can  be  parameterized  with  either  q  or  tjs  q.  Equations  (4—1)  and  (4-4)  are 
simple  functions  of  6  q  and  8 q.  The  converse  is  not  true;  the  6g  g  are  complex  func¬ 
tions  of  tjs  g.  The  6$  g  are  calculated  by  integrating  the  continuity  equation.  Each  S$  q  is 
the  sum  of  rj ^  g  for  different  segments.  The  sum  is  constrained  to  be  zero  at  the  end  of  the 
wire.  The  6  q  give  the  simplest  description  of  the  first  charge  basis  function. 
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The  higher  order  current  basis  functions  modify  the  current  in  the  interior  of  the 
segment  s  without  changing  the  current  at  the  ends  of  the  segment.  The  second  current 
basis  function  is  a  quadratic  polynomial  on  the  segment,  which  is  zero  at  the  ends  and  out¬ 
side  the  segment  (figure  4—1,  row  b).  The  third  basis  function  is  a  cubic  polynomial,  which 
is  zero  at  the  center  of  the  segment,  at  the  ends  of  the  segment,  and  outside  the  segment 
(figure  4—1,  row  c).  In  general,  the  n^  order  current  basis  function  is  an  n+1  order  poly¬ 
nomial,  which  is  zero  at  n— 1  interior  points  on  the  segment,  at  the  ends  of  the  segment, 
and  outside  the  segment.  However,  the  locations  of  the  n— 1  interior  zeros  is  not  obvious 
(figure  4—1,  row  d).  On  the  other  hand,  the  charge  basis  functions  (figure  4-1,  column  2) 
computed  from  the  continuity  equation,  are  proportional  to  the  Legendre  polynomials. 

In  the  general  case,  the  charge  basis  functions  are  Legendre  polynomials  on  the  inte¬ 
rior  of  the  segment  and  zero  outside  the  segment  (figure  4—1,  column  2).  The  charge  distri¬ 
bution  on  the  segment  is  approximated  by  N  charge  basis  functions: 

Qs(ts)  =  lfiOP0(ts)+Ni1>)snPA)' 

L  11  =  1 

Qs(‘s_j)  =  0  p0(Vj)- 

L 

The  current  distribution  is  calculated  from  the  continuity  equation  and  equa¬ 
tion  (4-5).  The  following  relations  are  identities  for  Legendre  polynomials  for  all  t: 

t 

PQ(p)dp  =  [PQ(t)  +  P1(t)]  for  n  =  0,  (4-6) 

-1 


(4-5) 
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Fn(p)dp = dr[p"+i(t)  -  Pn-i(t)1  for  n  >  0i 


(4-7) 


P  (t)  =  — —  —  (t2— l)n  for  n  >  0. 
n  2nn!  dtn 


(4-8) 


Equations  (4-6)  and  (4—7  for  n=l)  are  verified  by  direct  calculation.  In  the  general  case, 
equation  (4—7)  is  differentiated  and  simplified  with  the  equation  (4—8)  (Gradshteyn  and 
Ryzhik,  1980): 


1  1  jn+2 


2n  A(n— 1)!  dtn 


For  n  >  2,  this  equation  can  be  rewritten 


The  second  derivative  is  evaluated  to  give 


(4-9) 


(4—10) 


P"(t)  =  (2n+l)2^V+D^[t4(n+1)llt2+  (n+1)2(t2-1)-4(n+1)n](‘2-1)1 


(4-11) 
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The  three  terms  combine  to  give 


P  ft)  =  _ - _ - 

n  (2n+l)2n+1(n+l)!  dtn[ 


(n+l)(4n+2)(t2— l)n 


=  —  — (t2-l)n 
0n  ,  ,.nv  ' 

2  n!  dt 


(4-12) 


The  current  basis  functions  computed  from  the  continuity  equation  (3—3)  are 

Js  A  )  =  -  f^[pn+l  (‘4  )  -  Pn-A  )1  “a  for  ”  >  °> 

(4— 13  a) 
and 

J$  n(ty)  =  0  for  ft  s.  (4— 13b) 

These  current  basis  functions  modify  the  current  in  the  interior  of  the  segment  without 
changing  the  current  at  the  ends  of  the  segment. 

The  current  on  the  segments  s  and  s—1  is  approximated  by  N  current  basis  functions 


N-l 


and 


=  -‘A o;.-  1  f^WM-ViMvw 

n =  1  4n+^ 


Js^s-i^  -  ^1  +  ls-i^s  o  Vr 


(4— 14a) 


(4— 14b) 


Equation  (4—14)  models  the  current  on  the  wire  as  a  continuous  function  of  t.  The 
triangular  basis  function  creates  a  piecewise  linear  approximation  of  the  current.  The 
higher  order  basis  functions  represent  the  difference  between  the  current  distribution  and 
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the  piecewise  linear  approximation  of  the  current  distribution.  For  a  slowly  changing  func¬ 
tion,  the  approximation  error  goes  to  zero  as  the  segment  length,  L,  is  decreased. 

The  charge  basis  functions  model  the  charge  distribution  on  each  segment  as  a  con¬ 
tinuous  function.  However,  the  charge  distribution  at  the  junction  of  two  segments  may 
not  be  a  continuous  function.  In  nature,  this  problem  does  not  occur;  the  electric  field  at 
the  discontinuity  would  redistribute  the  charge  to  eliminate  the  discontinuity.  The  electric 
field  boundary  condition  should  create  a  continuous  charge  distribution.  The  charge  distri¬ 
bution  computed  for  the  dipole,  figures  2—8  and  2—9  in  section  2,  has  a  very  small  disconti¬ 
nuity  between  segments. 
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5  MULTIPOLE  POTENTIALS  OF  THE  BASIS  FUNCTIONS 


An  exact  model  for  the  continuous  charge  distribution  on  a  segment  requires  an 
infinite  sum  of  Legendre  polynomial  basis  functions.  In  our  numerical  model,  only  a  fixed 
number,  N,  of  charge  basis  functions  can  be  used  in  the  model.  An  error  in  modeling  the 
charge  distribution  introduces  errors  in  modeling  the  potentials.  The  error  in  the  potential 
is  discussed  for  three  different  regions:  the  far,  intermediate,  and  near  fields  of  the  segment. 
These  basis  functions  allow  the  errors  in  the  potentials  and  in  the  charge  distribution  to  be 
analyzed  and  modeled  in  detail. 

The  modeled  far,  intermediate,  and  near  fields  are  the  sum  of  the  fields  created  by 
the  individual  basis  functions.  Each  basis  function  contributes  differently  to  the  far,  inter¬ 
mediate,  and  near  fields  of  the  antenna.  The  contribution  to  the  three  regions  can  be 
examined  by  expanding  the  free-space  Green’s  function  in  a  1/r  expansion  (Jackson,  1975). 
In  this  section,  the  expansion  of  the  free-space  Green’s  function  is  limited  to  the  axis  of  the 
segment.  In  the  near  field,  the  basis  functions  are  multipole  expansions  of  the  charge  dis¬ 
tribution  on  a  segment.  In  particular,  the  basis  functions  represent  the  monopole,  dipole, 
quadrupole,  octupole,  etc.,  of  the  charge  distribution.  The  potentials  in  the  intermediate 
and  far  fields  depend  on  the  wavelength,  A,  of  the  exciting  frequency  and  segment  length. 
The  intermediate  and  far  field  terms  are  small  for  electrically  small  antennas.  This 
approach  is  compared  to  the  matrix  transformation  developed  by  Canning  (1992). 

Equation  (3 — 4),  without  the  constant  l/(47re0),  is  called  the  potential  integral: 


k  (r)  = 

s  nv  ’ 


P  (2z'/L> 
s  nv  '  '  i 


-jk  |  r— r' 


dV'. 


V' 


r— r' 


(5-1) 
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The  asymptotic  form  of  the  potential  integral  k  is  found  by  expanding  the  Green’s  func- 

s  n 

tion  in  a  1/r  series.  The  analysis  of  the  integral  is  simplified  by  limiting  its  evaluation  to 
the  z-axis  of  the  wire  and  setting  a  =  0.  The  potential  integral,  ks  n(z  )  on  the  z-axis,  is 
found  by  replacing  the  quantity  |  r-r'  |  in  equation  (5-1)  with  | z0— z  | ,  where  z^  is  the  dis¬ 
tance  to  the  center  of  the  charge  distribution,  Pn(2z/L),  on  the  segment  s  : 


fL/2 

Pn(2z/L) 

— L  /  2 


-jk(z0-z) 
-  dz. 

Z0-Z 


(5-2) 


If  |L/2|  <<  z0,  the  factors,  1/ 1 z0 — z  |  and  e^z,  in  the  free-space  Green’s  function  are 
expanded  in  a  power  series  of  z/z  and  jkz,  respectively: 


fL/2 

s  nM  = e_Jkz"  P„(2z/L) 

-L/2 


00  Jj  00 

k  =  0  zo  J  Lp  =  0 

for  z0  >>  L/2. 


kP2P 


dz 


(5-3) 


The  terms  with  like  power  of  zm  are  collected  together  as: 


K, 


s  n 


(z0) 


-jkz 


o  S 


m 

E 


7  P  k  P 


m-p+1 
m  =  0  Lp  =  0  p!z  o 


L/2 

Pn(2z/L)2md2, 

-L/2 


(5-4) 


where  m  =  p+k.  The  integral  can  be  simplified  by  making  a  change  of  variables  2z/L  =  t: 


K 


nW  = 


_  p-jLz 


o  E 

m  =  0 


;,7PkP(L/2)m+1 

“  -P+1 


p=0  p!z 


m- 


rl 


pn(t)t 


m 


dt. 


-1 


(5-5) 
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The  integral  factor  is  independent  of  frequency  and  segment  length.  The  integral  can  be 
evaluated  by  expanding  tm  in  terms  of  Legendre  polynomials: 


tm  =  s  (mp  (t); 

q  =  0 


(5-6) 


where  =  0  for  q  >  m.  The  integral  reduces  to 


1 


rl 


Pn(t)tmdt  =  s  c. 

«=» 


.m 


Pn(t)Pq(t)dt. 


-1 


(5-7) 


The  integration  of  the  orthogonal  Legendre  polynomials  gives 


Pn(t)tmdt  =  C^/(n  +  \). 


-1 


(5-8) 


For  m  <  n,  =  0,  the  above  integral  vanishes;  equatipn  (5—5)  simplifies  to 


Ks  n^z°^  =  e  jkZ°  S 


m  =n 


“  ;pkp(L/2) 

P.«  P!^“P+1 


m+T 


>IU  //  ,  1\ 

Cn/(“  +  2)- 


(5-9) 


At  very  large  distances  from  the  segment,  z^  >>  L,  most  of  the  terms  in  equation  (5—9)  are 
small.  For  a  fixed  p,  the  largest  term  occurs  for  m  =  n.  The  terms  with  m  >  n  involve 
higher  powers  of  l/z0.  The  sum  of  the  largest  terms  is 


!  \  jkzn  v  j PkP(L/2) ^n //  ,  1^ 

Ks  n(z»)  =  e  "  S  1 -  n-p+T  V<n  +  Z>  for 


|  11— D+  . 

p  =  0  P- z  0  v 


z0  >>  L. 


(5-10) 


3T 


Equation  (5—10)  exhibits  the  frequency  and  the  1/z  dependence  of  the  potential. 

The  term  with  p  =  0  is  independent  of  frequency  k;  it  is,  by  definition,  the  near  field 
part  potential  integral.  The  near  field  of  k  falls  off  as  [L/2z0]n"*"'*‘.  The  n^  basis  func- 
tion  is  the  n^  multipole  expansion  of  the  charge  distribution  on  a  segment.  The  intermed¬ 
iate  field  depends  on  the  frequency,  k,  but  falls  off  faster  than  l/z0,  n  >  p  >  0.  The  inter¬ 
mediate  field  represents  stored  energy;  they  do  not  radiate  energy  to  infinity.  The  term 
with  p  =  n  has  a  l/z0  dependence;  this  term  is  the  far  field  part  of  the  potential.  For  the 
far  field  term,  the  coefficient  of  the  (n/(n  +  i)  simplifies  to  j  nL/2[kL/2jn/(n!z0).  The 
factor,  kL/2,  simplifies  to  7tL/A.  These  results  are  summarized  in  figure  4—1  (Static  Part 
and  Far  Field  columns). 

Equation  (5—10)  has  a  simple  physical  interpretation.  For  a  static-charge  distribu¬ 
tion,  P  (t),  the  potential  cancels  at  large  distances  so  that  it  falls  off  as  l/z^~*~^  The  fre¬ 
quency  dependency  of  the  Green’s  function  prevents  the  potential  integral  from  canceling  in 
the  far  field  .  This  leads  to  frequency-dependent  terms  in  the  potential  integral. 

The  above  discussion  is  limited  to  the  radial  dependence  of  the  potential  integral. 

The  angular  dependence  of  the  near  field  (static  potential  integral)  is  Pn(cos<?),  where  9  is 

measured  relative  to  the  axis  of  the  segment.  For  the  case  n  =  0,  the  above  discussion 

needs  to  be  modified.  The  first  current  basis  function  yields  two  square-pulse  charge  basis 

functions,  which  form  a  dipole.  As  a  result,  the  potential  for  the  static-charge  distribution 

o 

created  by  the  first  current  basis  function  falls  off  as  l/z0. 


3.  If  the  segment  is  the  same  size  or  larger  than  the  wavelength,  the  opposite  effect  takes 
place.  In  this  case,  the  phase  term  causes  the  potentials  to  constructively  add  to  give  a  1/r 
term  for  some  of  the  basis  functions.  These  basis  functions  could  be  used  as  entire  basis 
functions;  however,  the  higher  order  basis  functions  may  not  fall  off  as  fast.  In  this  case, 
the  advantages  of  the  basis  functions  would  be  lost. 
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In  our  approach,  the  far  and  intermediate  fields  are  modeled  by  the  net  charge  dis¬ 
tribution,  PQ(ty),  on  each  segment.  The  high-order  basis  functions  model  the  near  field  as 
independent  degrees  of  freedom.  The  high-order  basis  functions  also  model  the  intermedi¬ 
ate  fields,  which  depend  on  the  factor,  kL/2  or  7rL/A,  in  equation  (5—10).  For  electrically 
small  antennas,  (h  <  A/4tt),  the  constant,  ttL/A,  is  small.  Therefore,  the  frequency- 
dependent  terms  in  the  potential  integral  are  small.  The  impedance  matrix  computed  from 
these  basis  function  has  small  frequency-dependent  terms.  The  properties  of  the  fields  for 
electrically  small  antennas  are  explicit  by  the  impedance  matrix. 

These  basis  functions  allow  the  error  to  be  discussed  in  detail.  There  are  two  types 
of  errors  in  this  model  of  the  charge-distribution  errors.  The  first  error  is  introduced  by 
using  a  finite  number  of  charge  basis  functions  to  model  the  charge  distribution.  The  sec¬ 
ond  error  is  caused  by  the  numerical  evaluation  of  the  charge  basis  functions  coefficient, 
77  .  Any  error  in  the  charge  distribution  creates  an  error  in  the  potentials.  However,  for 

these  basis  functions,  the  impact  of  the  errors  in  the  potentials  are  limited  to  the  near  field 

around  each  segment.  The  error  in  the  charge  basis  functions  coefficient,  77  ,  is  deter- 

s  n 

mined  by  the  electric  field  boundary  condition,  section  7.  The  detailed  discussion  of  the 
local  error  on  each  segment  is  given  in  section  8. 

These  basis  functions  create  a  structured  model  of  near,  intermediate,  and  far  fields; 
this  simplifies  the  impedance  matrix.  The  Canning  (1992)  used  a  matrix  transformation 
technique  to  add  structure  to  the  impedance  matrix.  The  matrix  transformation  creates 
oscillatory  basis  functions  that  model  the  charge  distribution  created  by  a  scattering  elec¬ 
tromagnetic  wave.  When  the  oscillations  in  the  charge  distribution  are  larger  than  a  wave¬ 
length,  the  oscillatory  basis  functions  have  directional  radiation  patterns.  The  directional 
radiation  pattern  simplifies  the  interaction  between  different  segments.  For  an  oscillation 
period  smaller  than  a  wavelength,  the  potentials  are  evanescent,  i.e.,  the  potentials  decay 
exponentially  in  the  radial  direction  from  the  segment.  The  matrix  transformation  will 
reduce  the  number  of  large  elements  in  the  impedance  matrix.  This  simplifies  the  solution 
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of  the  impedance  matrix.  However,  the  accuracy  of  the  solution  depends  on  the  accuracy 
of  the  algorithm  used  to  compute  the  initial  impedance  matrix. 

The  algorithm  developed  in  this  report  is  designed  to  accurately  compute  the  indi¬ 
vidual  elements  of  the  impedance  matrix.  The  basis  functions  and  boundary  condition 
combine  to  create  a  high-order  multipole  error  in  the  charge-distribution  model.  The 
admittance  at  resonance  and  antiresonance  depends  on  the  size  of  the  matrix  model.  The 
higher  order  basis  functions  simplify  the  impedance  matrix  and  improve  the  convergence  of 
the  algorithm. 

In  the  above  discussion,  the  wire  radius  was  not  included  in  the  1/z^  expansion  of 
the  potential.  However,  near  the  surface  of  the  segment,  the  finite  radius  must  be  included 
in  the  calculation  of  the  potential  integral.  The  numerical  calculation  of  the  potential  is 
discussed  in  section  6. 
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6  POTENTIAL  INTEGRAL 


The  convergence  of  a  numerical  algorithm  depends  on  the  accuracy  of  the  potential 
integral  calculation.  Approximations  in  the  potential  integral  calculation  can  cause  conver¬ 
gence  and  stability  problems  in  the  charge  distributions  at  the  end  of  the  wire.  These  prob¬ 
lems  are  discussed  by  Popovic  B.  D.,  M.B.  Dragovic,  and  A.  R.  Djordjevic  (1982).  The 

potential  integral,  k  on  the  surface  of  the  segment,  was  evaluated  by  Butler  (1975)  for 
s  n 

constant  and  linear  charge  distributions.  His  method  cannot  be  generalized  for  arbitrary 
basis  functions,  P  (s).  The  numerical  algorithm  developed  in  this  section  can  calculate  the 
potential  integral,  k  n,  to  any  desired  accuracy. 

In  section  5,  the  charge  distribution  was  on  a  line,  equation  (5—1).  In  this  section, 
the  charge  distribution  will  be  a  surface  charge,  l/(27ra)Pn(Lz'/2),  on  the  cylinder  of 
radius  a.  The  cylindrical  coordinate  system  (figure  6—1)  is  used  to  compute  the  potential 
integral.  The  axis  of  the  cylinder  is  the  z— axis  and  the  central  cross  section  of  the  cylinder 
is  in  the  x— y  plane.  The  potential  integral  is  evaluated  at  the  point  r=  ( p ,  tp,  z).  The 
kernel  has  rotational  symmetry  about  the  z-axis;  therefore,  the  equations  are  simplified  by 
setting  p=  0.  The  integration  is  over  the  surface  of  the  cylinder,  r'  =  (a,  tp' ,  z '): 


Kn(P>z)  =  — 

n  2tt 


■27T 


fL/2 


Pn(2z'/L)e  jklr  r'  I  dz'd tp' , 


-L/2 


T-T' 


where 


r—r'  |  =  [(z— z')2+/?2+a2— 2apcos^']'*'/2. 


(6-1) 


In  this  section,  the  subscript  s  will  be  dropped  from  the  k  (r)  and  the  charge  distribution, 

u  XI 

P  (t  ).  The  kernel  can  be  written  as  the  sum  of  two  parts:  the  singular  part  is 
n  s 
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The  evaluation  of  the  potential  integral  can  be  simplified  by  looking  at  the  kernel  of 
equation  (6—1).  The  value  of  the  potential  integral  depends  on  the  distance,  |  r—r'  | ,  from 
the  field  point  to  a  point  on  the  segment.  In  most  cases,  this  distance  from  the  segment  to 
the  field  point  is  large  compared  to  the  radius  of  the  wire.  The  integration  along  the  length 
of  the  segment  is  an  excellent  approximation  of  the  potential  integral,  «n;  the  angular  inte¬ 
gration  is  a  small  correction  to  the  potential  integral.  For  this  reason,  we  integrate  over 
the  length  of  the  segment  first,  then  use  a  Gauss-Chebyshev  integration  formula  to  inte¬ 
grate  around  the  circumference  of  the  segment.  A  sequence  of  Gauss-Chebyshev  integra¬ 
tion  formulas  is  used  to  compute  the  angular  integration.  Close  to  the  surface  of  the  seg¬ 
ment,  the  kernal  can  be  modified  to  improve  the  convergence.  The  angular  integration  can 
be  performed  to  any  desired  accuracy. 

Butler  (1975)  analytically  integrated  the  singular  term  in  the  kernel.  The  angular 
integration  was  expanded  in  power  series  and  integrated  term  by  term.  For  the  general 
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Legendre-polynomial  charge  distribution,  the  power  series  is  too  complicated  to  be  of  use. 
The  nonsingular  kernel  can  be  easily  integrated  with  a  sequence  of  Gauss-Chebyshev  inte¬ 
gration  formulas. 

This  numerical  method  is  illustrated  for  the  simple  case  n  =  0.  The  general  results 
are  given  in  appendix  A.  The  integration  along  the  surface  of  the  cylinder  yields 


•27T 

In 

1+r 

1— r 

Jo 

dr, 


(6-4) 


where  r  = - , 

R-L  =  /(z+L/2)2+p2— 2apcosy>'  +a2, 
and 

R2  = 

The  variable,  Rp  is  the  distance  from  the  bottom  edge  (a,  r  >  — L/ 2)  of  the  cylinder  to  the 
field  point.  The  variable,  R2,  is  the  distance  from  the  top  edge  (a,  r ,  L/2)  of  the  cylinder 
to  the  field  point.  The  argument  of  the  natural  logarithm  in  equation  (6-4)  can  be  calcu¬ 
lated  very  accurately  by  using  one  of  several  different  algebraically  equivalent  forms,  equa¬ 
tions  (A— 13),  (A— 14),  and  (A— 15)  in  appendix  A.  The  potential  near  the  surface  of  the 
wire,  equation  (A— 15),  can  be  simplified  to  improve  the  numerical  convergence.  The  term 
ln(p2— 2a/3COsr+a2)  in  equation  (A— 15)  is  the  potential  of  an  infinite  line  charge  at  <p'  and 
p  =  a.  The  angular  integral  of  ln(p2— 2apcosr +a2)  gives  the  potential  for  an  infinite  cylin¬ 
der  ln(p2)  for  p  >  a  and  ln(a2)  for  p  <  a.  The  details  are  given  in  appendix  A. 

The  singular  part  of  equation  (6—4)  can  be  integrated  to  give  exact  results,  equa¬ 
tions  (6— 5a)  and  (6— 5b).  Equations  (6-5a)  and  ( 6-5b )  are  not  approximations. 
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-j  (z+L/2+R  )(R  — z+L/2) 

Kft(p,z)  =  —  In  - 1 - - -  dtp'  for  p  >  a  and  |z|  <  L/2, 

2tJo  I-  p2  J 


(6— 5  a) 


r2  7T 

«n(p,z)  =  ~  ln 

2dn 


(z+L/2+R  )(R  -z+L/2)' 

- - - - -  dip'  for  p  <  a  and  jz|  <  L/2 

a2 


(6— 5b) 


The  potential  at  the  end  of  the  segment  |  z|  =  L/2  is  discussed  in  appendix  A. 

The  integral  around  the  cylinder  can  be  simplified  by  making  the  substitution 
COS^'  =  x. 


«q(p,z) 


i 

7 r  1— r(x 


(6-6) 


This  integral  can  be  evaluated  by  using  the  Gauss-Chebyshev  integration  formula.  The 
Chebyshev  polynomials  are 


T»t(x)  =  cos(N  cos  i(x))J 


where  — 1  <  x  <  1.  The  quadrature  points  are  the  roots  of  T^(x): 


N  cos  ^(x.)  =  +  i7r  for  0  <  i<  N— 1 


x.  =  cos[7r(i+l/2)/N]. 


45 


In  the  Gauss-Chebyshev  integration  formula  each  quadrature  point  has  a  weight 


w-  =  7t/N. 

The  Gauss-Chebyshev  integration  formula  for  our  integral  is 

r1 

f(x) 

-1 

The  constant  l/7r  in  equation  (6—6)  can  be  expressed  in  terms  of  the  total  charge  on  the 
wire,  Q  =  L. 


dx 


y  1- y, 


N-l 
=  E 


Wif(xj). 


i  =0 


=  E_1 

LN  i=0 


In 


'1  +  T(Xj)' 


(6-7) 


The  angular  integration  in  equation  (6—7)  has  a  simple  interpretation.  The  con¬ 
stant  surface  charge  density  is  approximated  with  2N  line-charge  distributions  with  a  con¬ 
stant  charge  per  unit  length.  The  total  charge  on  each  line  is  Q/(2N).  The  line  charge 
distributions  are  located  at  =  7r(i+l/2)/N  and  p'  =  a;  the  symmetry  in  the  angular 
distribution  allows  the  number  of  terms  in  the  sum  to  be  reduced  to  N.  Equation  (6—7)  is 
equivalent  to 


\(p,z)  = 


L/2 


-0-  E  1 


-L/2  ^  *  ==0  [(z— z,)2+/?2+a2-2apcos^.] 


1/2 


dz. 


(6-8) 


When  N  =  1,  equation  (6—8)  is  the  reduced  kernel  approximation.  In  the  reduced 
kernel  approximation,  |  r—  r'  \  is  replaced  by  [(z-2/)2+p2+a2j.  Equation  {6—8)  is  a  generali- 
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zation  of  the  reduced  kernel  approximation.  The  angular  integration  of  the  potential  inte¬ 
gral  is  the  same  for  the  general-charge  distributions,  Pn(2z'/L). 

The  error  in  the  angular  integration  can  be  calculated  by  using  the  symmetry  of  the 
2N  line  charges  and  a  spherical  harmonic  expansion  (Jackson,  1975)  of  k °(p,z).  The  spheri¬ 
cal  harmonic  expansion  is 

S  S  Bh:-(,+1)Yh(W),  (6-9) 

i=0  m=-l 

which  is  valid  in  a  charge-free  region,  r  >  ((L/2)2  +  a2)^.  The  variables  r,  0,  and  <j> 
define  a  spherical  coordinate  system;  r  is  the  distance  from  the  origin,  6  is  the  angle  from 
the  z-axis,  and  <(>  is  the  angle  measured  around  the  z-axis  =  0  in  x— z  plain).  The  m  —  0 
terms  are  cylindrically  symmetric.  The  mi  0  terms  have  an  m— fold  rotational  symmetry 
and  vanish  on  the  z-axis. 

The  rotational  symmetry  of  the  kernel  ensures  that  all  the  Gauss-Chebyshev  Integra- 

g 

tion  formulas  are  exact  on  the  z-axis.  The  potential  integral,  n  on  the  z-axis,  uniquely 
determines  the  value  of  all  the  cylindrically  symmetric  terms;  they  are  the  same  for  all 
Gauss-Chebyshev  integration  formulas.  The  nonsymmetric  terms  represent  the  error  in  the 
Gauss-Chebyshev  integration  formulas. 

The  2N  line-charge  distributions  and  the  spherical  harmonic  expansion  must  have 
the  same  rotation  symmetry.  The  2N  line-charge  distributions,  Pn(2z'/L),  are  symmetric 
under  rotations  of  27t/2N  radians.  Consequently,  the  leading  error  term  has  m  =  2N, 
where  l  >  m.  The  value  of  l  depends  on  the  symmetry  of  the  charge  distribution  on  the 
segment.  The  charge  distribution  Pn(0)  is  even  for  even  n;  the  potential,  is  even 

if  l  is  even  and  l>  2N.  The  charge  distribution,  Pn(0),  is  odd  for  odd  n;  the  potential, 
Y im (0,0),  is  odd  if  l  is  odd  and  l  >  2N+1.  The  radial  dependence  of  the  first  error  term  is 
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r— (2N+1)  £oj.  gven  n  an(|  r— (2N+2)  £oj.  0(^  n  accuracy  0f  the  angular  integration 

increases  quickly  with  the  order  N  of  the  Chebyshev  quadrature  formula. 

This  error  analysis  must  be  modified  if  r  <  ((L/2)2  +  a2)"^.  In  this  case,  the 
potential  is  expanded  using  modified  Bessel  functions  (Jackson,  1975).  In  a  cylindrical 
coordinate  system: 


= 


S  (A  (k)cos(kz)  +  B  (k)sin(kz))K  (kp)e  Mlk. 
v=0  v 


(6-10) 


This  expression  is  valid  for  p  >  a.  The  variable  k  is  a  continuous  parameter;  k  is  not  the 
wave  number.  The  v  =  0  terms  are  cylindrically  symmetric.  The  error  term  has  a  2N 
rotational  symmetry,  v  =  2N,  4N,  6N,  etc.  In  this  case,  the  error  is  a  continuous  function 
of  k.  The  radial  dependence  is  determined  by  the  K^(k p)  factor.  The  factor,  K^(k p),  has 
the  limits: 


2 

Jp. 


i  V 


for  v  ±  0  and  kp  <<  1, 


(6— 11a) 


Kz/k/5)  "  |_2kpJ  e 


7 T 


J/2  -k„r 


1  +  0 


'  1 
Tp] 


for  k p  »  1. 


(6-1  lb) 


Equation  (6-lla)  applies  near  the  surface,  p  ~  a  and  k  <<  1/a.  Equation  (6— lib)  applies 
at  larger  distances,  p  ~  L/2,  where  k  >>  2/L.  For  most  values  of  k,  the  error  terms  fall  off 
as  p  near  the  surface  of  the  wire. 

The  angular  integration  can  be  computed  with  a  sequence  of  Chebyshev  quadrature 
formulas  with  N  =  1,  3,  9,  27,  81.  This  sequence  requires  the  minimum  number  of  new 


48 


quadrature  points^  for  each  new  iteration.  The  error  in  the  iterations  N  =  1,  3,  9,  27,  81 
falls  off  as  1/r  (or  1/p)  to  the  power  3,  7,  19,  55,  163,  respectively.  The  power  of  the  r 
more  than  doubles  for  the  second  and  third  iteration;  it  almost  triples  for  the  fourth  and 
fifth  iteration. 

The  error  in  each  integration  can  be  estimated  by  using  the  previous  iteration.  The 

s 

normalized  error^  in  the  p  iteration  of  the  static  potential  integral,  «n(r,#,0)p,  is 


exact 


(6-12) 


g 

where  /in(ri^0)exaC|j  *s  the  exact -potential  integral.  The  error  in  the  sequence  of  iterations 
obeys  the  relationship: 


2 

error  <  error^.  (6—13) 

The  normalized  iteration  difference  (NID)  is 


NID  (p-l,p)  = 


Kn(r^>0)p-l~Kn(r’^°)i 


(6-14) 


The  NID  can  be  expressed  in  term  of  the  error^^ 


NID  (p— l,p)  ~  error  .—  error  . 

p— 1  p 


(6—15) 


The  second  term,  erro^,  in  the  above  equation  is  much  smaller  than  the  first  term;  it  can 


4.  L.  Koyama  pointed  out  that  the  Gauss-Chebyshev  integration  formulas  have  some  com¬ 
mon  evaluation  points. 
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In  terms  of  the 


be  neglected.  The  error  in  the  p^ 
NID(p— l,p),  the  error^  is 


iteration  is  less  than  error 


2 

p-r 


erroip  <  NID(p-l.p)^. 


(6-16) 


The  error  in  the  present  iteration  is  the  less  than  square  of  the  difference  between  the  present 

—8 

and  previous  iterations.  The  iteration  sequence  can  be  stopped  when  the  NID  is  10  for 
double  precision  and  10  ^  for  single  precision. 

The  above  conclusions  are  verified  with  a  computer  calculation  of  the  static 
potential  integral  at  a  distance,  L,  from  the  center  of  a  segment  with  a  length,  L,  and 
radius,  1/2.  The  NID  for  n  =  0  and  1  was  calculated  for  the  range  0  <  0  <  7r/2;  the  largest 
value  for  this  range  is  given  in  table  6—1.  For  n  >  1,  the  NID  is  normalized  with  the 
potential  on  the  z-axis  is 


NIDn(p-l,p) 


for  n  >  1. 


(6—17) 


The  largest  NID  for  n  >  1  is  given  in  table  6—1.  The  NID  in  the  second  row  is  less  than 
the  square  of  NID  in  the  first  rows.  The  NID  in  the  third  row  is  limited  by  the  precision  of 
the  computer. 

In  section  5,  we  discussed  the  functional  form  of  the  potential  at  large  distances 
from  the  segment.  This  algorithm  allows  the  static  potential  integral,  ks  ,  to  be  computed 

o  XI 

o 

near  the  surface  of  a  wire.  Figure  6— 2a  plots  the  static  potential  integral,  ,  along  the 
surface  of  the  source  segment  and  the  nearest-neighbor  segments.  The  source  segment  is 
centered  at  x  =  0,  with  a  length  of  10  and  a  radius  of  1/2.  Figure  6 -2b  uses  an  expanded 

g 

scale  for  the  static  potential  integrals,  k  n,  on  the  surface  of  the  next  two  segments.  This 
plot  shows  that  the  multipole  static  potential  integrals,  k  ,  falls  off  very  fast. 

5  XI 
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Figure  6-2.  Potential  integral  on  the  surface  of  a  wire  produced  by  a  source  on  segment  l, 
where  (a)  is  the  potential  near  source  and  (b)  is  the  potential  on  the  second  and  third  segments 
(length,  L=100,  and  radius,  a=1/2). 
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Table  6—1.  Normalized  iteration  difference  in  the  angular  integration  on  a 
spherical  shell  with  the  radius  of  10  for  (a)  and  100  for  (b). 


Iteration 

Monopole 

Dipole 

Quadrupole 

Octupole 

Difference 

1-2 

2x10  3 

3x10  3 

3x10  3 

4x10  3 

- 

2-3 

7xl0-9 

3xl0-8 

8xl0-8 

2xl0“7 

3 -A 

r  ,n— 16  K  in— 12  o  in— 12  ir>— 12 

5x10  5x10  3x10  6x10 

(a)  Segment  length  is  10  and  radius  is  1/2. 

Iteration 

Monopole 

Dipole 

Quadrupole 

Octupole 

Difference 

1-2 

2x10  5 

4x10  5 

6x10  5 

lxlO-4 

2-3 

7xl0~15 

2xl0~n 

r— 1 

1 

o 

t— H 

X 

00 

2xl0— 13 

3-4 

5xl0-15 

,  in— 15 
5x10 

3x10  16 

6xl0-17 

(b)  Segment  length  is  100  and  radius  is  1/2. 


Figures  6— 3a  and  6— 3b  show  the  potentials  for  a  higher  aspect-ratio  segment  with  a 

length  of  100  and  radius  of  1/2.  The  aspect  ratio  is  the  segment  length  divided  by  the 

segment  radius.  Figure  6— 3a  shows  the  static  potential  integral,  ks  ,  for  the  source  seg- 

s  n 

o 

ment  and  the  two  adjacent  segments.  Figure  6— 3b  shows  the  static  potential  integral,  «  , 

s  n 

of  the  next  two  segments.  The  plots  show  that  the  high-order  static-potential  integral, 

g 

k  ,  is  large  only  on  segments  near  the  source  segment. 

g 

The  static  potential  integral,  k  ,  for  the  two  different  aspect -ratio  segments  have 

e>  Ii 

the  same  shape  as  the  Legendre  polynomial  charge  distribution.  Therefore,  the  potentials 
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POTENTIAL  INTEGRAL 


monopole 


dipole 


quadbapole 


OCTOPOLE 


D,  stance  along  subfaceofw.be 


monopole  I 
DIPOLE 
QUADRAPOLE 
OCTOPOLE 


l^wwwWggiSBB^ — w 


distance  along  surface  of  wire 

(b) 


w  t  c 


should  have  the  same  values  at  the  end  of  the  segment  where  P  (1)  =  1.  The  static 

s 

potential  integrals,  /c  ,  for  the  high  aspect-ratio  segment  are  nearly  the  same  at  the  end  of 

o  XI 

the  segment.  For  the  short  segment,  the  difference  in  the  static  potential  integrals,  ks  ,  is 

s  n 

large  at  the  end  of  the  segment.  In  addition,  the  slope  of  the  potentials  at  the  end  of  the 
segment  is  much  higher  for  the  high  aspect-ratio  wire.  The  shape  of  the  potentials  at  the 
end  of  the  wire  is  determined  by  the  details  of  the  local  charge  distribution  at  the  end  of 
the  segment. 

The  charge  distribution  on  the  short  segment  changes  very  fast  compared  to  the 

wire  radius.  Any  details  in  the  charge  distribution  smaller  than  a  wire  radius  are  averaged 

out  by  the  potential  integral.  For  high  aspect-ratio  segments,  the  charge  distribution 

changes  very  slowly.  For  high  aspect-ratio  segments,  the  values  of  the  potential  integral 

will  be  almost  the  same  near  the  end  of  the  segment  (t  =1).  The  sharp  falloff  in  the 

s 

potential  at  the  end  of  the  segment  is  caused  by  discontinuity  in  the  charge  distribution  at 
the  end  of  the  segment.  The  sharp  drop  is  limited  to  a  distance  of  about  six  wire  radii. 
This  sharp  drop  in  the  potential  can  affect  the  accuracy  of  the  boundary-condition  calcula¬ 
tion.  Details  smaller  than  a  wire  radius  cannot  be  resolved  by  this  approach. 

c  C 

The  shape  of  the  vector  potential  integral,  n_p  also  depends  on  the 

shape  of  the  current  basis  functions.  In  contrast  to  the  charge  basis  functions,  the  current 


basis  function  is  zero  at  the  end  of  each  segment.  This  eliminates  the  sharp  falloff  in  the 

s  s 

shape  of  the  vector  potential  integral,  n_p  Figures  6-4  and  6—5  show  the 

s  s 

shape  of  the  vector  potential  integral,  n  ,  j—  k$  n_p  for  second  and  third  current  basis 
function.  Figure  6-4  is  for  a  segment  with  a  length  of  10  and  a  radius  of  1/2.  Figure  6—5 


is  for  a  longer  segment  with  a  length  of  L  =  100  and  a  radius  of  1/2. 


The  plot  shows  that  the  higher  multipole  basis  functions  contribute  very  little  to  the 

g 

scalar  potential  integral,  k  ,  on  nonadjacent  segments.  The  falloff  in  the  potential  inte- 

o  XI 

s 

gral,  k  ,  determines  the  magnitude  of  the  off-diagonal  terms  in  the  impedance  matrix. 

o  II 
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A  detailed  discussion  of  the  individual  terms  and  error  in  the  impedance  matrix  is  given  in 
section  8. 


PART  B:  FREQUENCY-DEPENDENT  POTENTIAL  INTEGRAL 

The  frequency  part  of  the  kernel  can  be  numerically  integrated  with  the  same  tech¬ 
nique  used  in  part  A.  In  this  case,  the  integral  of  each  line-charge  distribution  is  evaluated 
by  Gauss-Legendre  integration.  The  annular  integral  will  be  integrated  by  using  the  same 
sequence  of  Gauss-Chebyshev  integration  formulas.  The  number  of  Gauss-Legendre  quad¬ 
rature  points  needed  in  the  integration  of  the  line-charge  distributions  depends  on  the  num¬ 
ber  of  segments  per  wavelength.  In  this  report,  we  will  assume  two  or  more  segments  per 
wavelength.  The  error  analysis  used  in  part  A  does  not  apply  to  the  integration  of  this 
part  of  the  kernel. 

The  integration  of  the  frequency-dependent  part  of  Green’s  function  gives  a  complex 
result.  The  real  and  imaginary  part  of  the  integral  are  considered  separately.  The  imagi¬ 
nary  part  has  the  form: 


(6-18) 

The  integrand  and  its  derivatives  are  continuous  for  all  |  r—r'  [ .  This  integral  can  be  accu¬ 
rately  evaluated  with  a  Gauss-Legendre  integration  formula. 

For  electrically  small  antennas,  the  integral  will  have  large  numerical  errors  for  the 
high-order  basis  functions.  An  error  in  the  evaluation  of  this  integral  will  introduce  an 


=  —  k 


N— 1  rL/2 
£ 

LN  i  =0 


Pn(2z '  /L)sin(k  |  r—r'  \  )/(k  |  r—r'  |  )dz ' . 
-L  /  2 
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error  in  the  feed  point  resistance  of  an  electrically  small  antenna.  Equation  (6—18)  can  be 
optimized  for  low  frequencies  by  using  an  alternative  form: 


n  LN  i  =0 


,L/2 

^  sin(k|  r-r’ | )" 

-L^  + 

Pn(2Z7L) 

dz' 

(kk-r’|). 

- 

-L/2 

- 

(6-19) 


The  constant  terms,  L<5^j,  cancel  for  the  two  square-charge  pulses  created  by  the  triangular 
current  basis  function.  This  cancellation  would  give  a  large  roundoff  error  for  electrically 
small  antennas.  The  impedance  matrix  is  computed  with  the  right  side  of  equation  (6—19). 
The  constant  term,  hfl,  is  computed  as  a  separate  term;  the  l/|  term  is  added  to  the 
potentials  used  to  compute  the  Js  q  and  Jg  ^  terms  in  the  impedance  matrix.  The  differ¬ 
ence  between  the  kernel  and  unity  must  be  computed  with  an  infinite  series  to  achieve  the 
desired  accuracy. 

The  real  part  of  the  frequency  potential  integral,  has  the  form: 


Re^(r)  =  - 


N-l  rL/2 
£ 


LN  i=0 


Pn(2z'/L)sin:‘ 


k  I  r-r' 


/(I  T~T'  I )dz' . 


-L/2 


(6-20) 


The  integrand  is  continuous,  but  its  derivative  does  not  exist  at  [  r-r'  |  =  0.  The  integra¬ 
tion  of  equation  (6—20)  is  simple  for  large  |  r—r'  |  >  L.  For  small  |  r—r'  \  <  L,  the  integra¬ 
tion  interval  is  divided  into  two  smooth  parts.  The  dividing  line  is  the  point  where  |  r—r'  \ 
is  minimum.  On  the  source  segment,  a  7-place  accuracy  can  be  achieved  by  using  a  42- 
point  Gauss-Legendre  integration  formula  for  each  interval.  On  the  other  segment,  a  14- 
place  accuracy  can  be  achieved  with  42-point  Gauss-Legendre  integration  formula. 
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7  BOUNDARY  CONDITION 


The  boundary  condition  for  a  perfect  conductor  requires  the  electric  field  to  be  per¬ 
pendicular  to  the  surface  of  the  conductor.  Therefore,  the  parallel  component  of  the  elec¬ 
tric  field,  E||,  must  vanish  on  the  surface  of  the  conductor.  For  thin-wire  antennas,  the 
boundary  condition  is  only  imposed  on  the  electric  fields  parallel  to  the  length  of  the  wire. 
No  boundary  condition  is  imposed  upon  the  electric  fields  perpendicular  to  the  length  of  the 
wire.  This  approximation  introduces  an  error  in  the  angular  charge  distribution  on  the 
wire.  Our  discussion  is  limited  to  antennas  excited  at  a  single  feed  point. 

In  the  method  of  moments,  the  boundary  condition  is  an  approximation  of  the  ideal 
boundary  condition.  The  approximation  of  the  ideal  boundary  condition  introduces  a  mod¬ 
eling  error  in  the  boundary  condition.  The  boundary  condition  is  calculated  from  the  sca¬ 
lar  and  vector  potentials.  The  scalar  and  vector  potentials  are  implicitly  a  function  of  all 
the  charge  and  current  distributions  on  all  the  segments.  The  modeling  error  in  the  charge 
distribution  will  induce  an  error  in  the  boundary  condition.  Likewise,  the  modeling  error  in 
the  boundary  condition  also  induces  an  error  in  the  charge  distribution.  In  addition,  the 
numerical  evaluation  of  the  boundary  conditions  will  introduce  errors  in  the  boundary  con¬ 
dition.  The  interaction  of  these  errors  on  the  single  segment  is  discussed  in  this  section. 
The  general  discussion  is  given  in  the  next  section. 

The  boundary  condition  developed  in  this  section  is  a  generalization  of  the  pulse¬ 
testing  boundary  condition.  The  pulse-testing  boundary  condition  integrates  the  product  of 
E||  and  a  testing  function.  The  boundary  condition  is  imposed  by  requiring  the  tested  elec¬ 
tric  field  to  vanish.  The  testing  function,  w(s),  is  typically  a  triangular  or  square-testing 
function.  In  this  report,  we  only  consider  square-testing  functions.  The  square-testing 
function,  w(z( t^)),  is  1  between  the  points  z{ t^)  and  z( tg);  it  vanishes  for  all  other  points: 
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<u/,Ei|>  = 


(7-1) 


r*(  t  ) 

2 

w(z )  E(z)-dz. 

The  boundary  condition  is 
<u/,E||>  =  0. 

In  our  case,  N  square-testing  functions  for  each  segment  are  required  for  the  N  charge  basis 
functions. 

For  square-testing  functions,  the  right  side  of  equation  (7—1)  is  the  integral  of  Ej| 
between  the  points  z{X^)  and  z^).  This  integral  is  the  electromotive  force  (EMF)  between 
the  points  t-^  and  The  EMF(tpt2)  is  the  energy  gained  by  moving  a  unit  charge  from 
point  tj  to  point  t2  on  the  surface  of  the  wire.  Ideally,  the  energy  gained  by  moving  a  unit 
charge  on  the  surface  should  be  zero  at  all  points.  The  EMF  can  be  expressed  in  terms  of 
the  potentials: 


EMF(y2)  =  ^(tp  —  <j>(t2)-iw 


A(z)  •  d  z . 


4* ,) 


(7-2) 


The  width  of  the  testing  function  does  not  need  to  be  the  same  for  each  testing  func¬ 
tion.  Michalski  and  Butler  (1983)  used  variable-width  testing  functions  to  analyze  scatter¬ 
ing  from  a  finite-length  cylinder.  The  boundary  condition  for  N  basis  functions  could  be 
imposed  by  using  N  testing  functions  for  each  segment.  The  functional  form  of  the  bound¬ 
ary  condition  error  on  each  segment  is  controlled  by  the  width  of  the  N  testing  pulses. 
However,  the  magnitude  of  the  numerical  error  on  the  segment  is  unknown.  The  functional 
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form  and  the  magnitude  of  the  error  can  be  derived  by  interpreting  the  boundary-condition 
equations  in  a  different  way. 

Instead  of  looking  at  the  EMF  for  individual  testing  pulse,  the  EMF  should  be  inter¬ 
preted  as  a  continuous  function  on  the  entire  antenna  being  modeled.  The  EMF(ty)  on  the 
segment  is  the  line  integral  from  the  center  of  the  feed  point  to  the  point  ty  on  the 
segment  f: 


EMF(ty)  = 


*(‘  /) 

E(z)  -dz  =  4>(0)  —  <j)(ty)  —joj 


z(  t  y  ) 
A(z)-dz . 

0 


(7-3) 


The  integration  of  the  vector  potential  is  a  line  integral  along  the  surface  of  the  antenna 
wires.  Ideally,  the  EMF(ty)  on  the  antenna  should  be  a  constant  equal  to  half  of  the  input 
voltage. 

The  EMF(ty)  function  on  a  segment  is  approximate  by  a  sum  of  Legendre  polyno¬ 
mials. 


N— 1 

EMF(t,)=  E 
•'  k=0 


0/kPt,l/(' 


(7-4) 


The  N  boundary  conditions  for  each  segment  are  imposed  by  setting  the  expansion  con¬ 
stants  equal  to 


“/  0  “  2Vin 


a/k  =  ®  for  0  <  k  <  N. 


(7-5) 


The  first  error  term  in  this  boundary  condition  is  PN(S)- 
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The  expansion  coefficient,  a^}  is  evaluated  by  computing  the  integral 


°/k  =  (k+Z> 


P  k(t/)EMF(t^)dt^- 


(7-6) 


This  integral  can  be  evaluated  very  accurately  with  an  M-term  Gauss-Legendre  integration 
formula,  where  M  >  N: 


where 


M— 1 

aA= ,  ^kiEMF(ti/)> 

4  i  =  (k+j)wiPk(‘i/)- 


(7-7) 


The  points  t-  are  the  quadrature  points  and  w^  are  the  quadrature  weights.  If  M  >  N, 
equation  (7—7)  gives  extra  expansion  coefficients.  These  extra  expansion  coefficients, 
to  represent  the  modeling  error  in  the  boundary  condition. 

The  boundary-condition  error  is  caused  by  using  a  finite  number  of  basis  functions 
to  model  the  charge  distribution.  The  charge  distribution  on  each  segment  is  modeled  by 
the  charge  basis  functions,  Pq^)  through  P^T_^(t j).  The  first  error  term  in  the  charge 
distribution  is  P^(t j)  charge  basis  function  contributes  an  error,  /t^(t ^),  to 

the  scalar  potential  in  the  boundary  condition.  The  «^(t j)  error  has  the  same  shape  as 
P^(ty)  on  the  source  segment.  The  high-order  charge  basis  functions,  P^(ty),  P N+l(t/)> 
etc.,  to  first  approximation  will  cause  errors,  pN(t/)>  Pn+i^/)’  etc'’  *n  tbe  boundary 
condition.  These  errors  could  affect  the  numerical  evaluation  of  the  boundary  condition. 

The  numerical  evaluation  of  the  integral  introduces  some  error  in  the  boundary 
condition.  The  M-point  Gauss-Legendre  integration  formula  is  exact  for  polynomials  with 


62 


a  degree  less  than  2M.  Thus,  the  first  error  term  in  the  evaluation  of  the  expansion  coeffi- 
dent,  ag,  is  caused  by  t  or  P ^(ty)  errors  in  the  EMF(ty)  function.  The  modeling 
error  terms  P^t  y)  to  P£M-l^f^  271  EMFfij)  function  have  no  impact  on  the  numeri¬ 
cal  calculation  of  the  coefficient  a.g.  As  a  consequence,  the  error  in  the  far  field  is  very 
small.  In  general,  the  expansion  coefficient,  a^,  is  exact  for  polynomials  with  a  degree  less 
than  2M— k.  The  error  in  the  expansion  coefficient,  is  larger;  however,  this  error  is  a 
high-order  multipole  moment.  The  high-order  multipole  moments  are  large  only  in  the 
near  fields.  The  error  terms  in  the  numerical  boundary  condition,  a-^  through  ,  are  a 
good  estimate  of  the  error  in  the  boundary  condition. 

The  above  summary  of  the  error  analysis  assumes  that  the  potential  integral  is  a 
smooth  function  that  can  be  modeled  with  low-degree  polynomials.  This  is  true  for  the 
interior  of  the  segment;  however,  it  is  not  true  at  the  end  of  the  segment.  The  distance 
from  the  peak  in  potential  integral  ks  ^  to  the  end  of  the  segment  is  6.7  wire  radius.  This 
sharp  drop  would  introduce  a  large  error  in  the  evaluation  of  the  expansion  coefficient, 
0:^^.  The  Gauss-Legendre  integration  with  M  =  4  requires  a  segment  with  an  aspect  ratio 
larger  than  85.  Table  7—1  gives  the  minimum  aspect  ratio  for  different  order  Gauss- 
Legendre  integration  formulas.  The  minimum  aspect  ratios  given  in  table  7—1  avoid  the 
sharp  drop  near  the  end  of  the  segment;  the  first  Gauss-Legendre  point  is  a  distance  of  6 
wire  radius  from  the  end  of  the  segment.  The  aspect  ratio  of  the  segment  limits  the  num¬ 
ber  basis  allowed  on  each  segment. 

The  end  effect  could  cause  an  error  in  the  electric-field  boundary  condition  at  the 
junction  of  two  segments.  However,  if  the  charge  distribution  is  continuous,  the  end  effects 


5.  The  general  form  for  the  error  term  in  the  integrand  of  the  M-point  quadrature  formula 
is  snPjyj-(s).  This  error  term  is  zero  at  the  quadrature  points.  The  factor  s11  can  be 

expressed  as  the  sum  of  Pq(s)  to  P-^(s).  The  first  term  to  have  a  non-zero  integral  is 

2  2M 

PM(s);  therefore,  the  first  error  term  in  the  integral  has  the  power  s 
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Table  7—1.  Minimum  aspect  ratio. 


Order 

Gauss-Legendre 

Integration 

Aspect  Ratio 

4 

85 

5 

130 

6 

180 

7 

240 

8 

300 

9 

375 

cancel  out.  At  the  free  end  of  a  wire,  the  end  effect  will  cause  an  error  in  the  boundary 
condition;  this  error  is  limited  to  a  6  radius  distance,  In  addition,  at  the  free  end  of  a  wire, 
several  high-order  polynomials  are  required  to  model  the  charge  distribution.  These  high- 
order  polynomials  may  introduce  errors  in  the  evaluation  of  the  boundary  condition.  The 
numerical  results  in  section  2  show  that  the  boundary  condition  is  very  stable  at  the  end  of 
the  wire. 

This  boundary  condition  is  simple  in  antennas  without  closed  loops.  In  a  loop,  the 
point  t j  can  be  reached  by  using  two  different  integration  paths.  The  EMF(t j)  needs  to 
be  independent  of  the  integration  path.  The  scalar  part  of  the  EMF  is  independent  of  the 
integration  path;  however,  the  integral  of  the  A  vector  depends  on  the  integration  path. 
The  integral  of  the  A  vector  on  a  closed  loop  is  evaluated  with  Stokes’  theorem, 

o  A(z)-dz=  (p,  (7-8) 


where  <j>  is  the  total  magnetic  flux  passing  through  the  loop.  In  a  closed  perfect  conducting 
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loop,  the  magnetic  flux  through  the  loop  is  zero.  The  EMF  boundary  condition  requires 
the  total  magnetic  flux  through  each  loop  to  vanish. 

The  above  arguments  indicate  that  the  coefficients,  a  ,  are  determined  very  accu¬ 
rately  for  a  smooth  charge  distribution.  Near  the  end  of  a  wire,  the  error  should  be  a  high- 
order  multipole  moment.  The  magnitude  of  the  error  in  the  EMF  boundary  condition  is 
easily  calculated. 
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8  IMPEDANCE  MATRIX 


In  this  section,  the  basis  functions  and  boundary  condition  are  used  to  construct  the 
impedance  matrix  for  the  method  of  moments.  The  boundary-condition  errors  can  be  esti¬ 
mated  by  computing  extra  terms  in  the  impedance  matrix.  The  errors  and  approximations 
used  in  constructing  the  impedance  matrix  are  discussed  in  detail.  For  the  dipole  problem, 
the  rows  and  columns  can  be  arranged  to  eliminate  the  need  for  pivoting  in  the  solution  of 
the  impedance  matrix.  The  Legendre-polynomial  basis  functions  are  compared  to  the 
power  series  basis  functions  used  by  Djordjevic  et  al. ,  (1990).  A  matrix  transformation  is 
developed  to  compare  the  EMF  boundary  conditions  to  the  pulse-testing  boundary  condi¬ 
tions.  The  EMF  boundary  condition  produces  a  matrix  with  smaller  off-diagonal  terms. 

The  impedance  matrix  relates  the  m^  boundary  condition  expansion  coefficient 
th 

a/ m  on  segments  /to  the  n  basis  function  coefficient  xs  n  on  segments  s: 


OLr  =  S  E  Z  r  y 
/m  fs  mnAs  n 

J  s  n  J 


(8-1) 


For  n  >  0,  yg  n  is  the  coefficient,  7?g  n,  for  the  charge  basis  functions  on  each  segment  s. 
For  n  =  0,  yg  g  is  the  coefficient,  5$  q,  of  the  current  basis  function.  Note,  the  free  end  of 
the  wire  does  not  require  a  current  basis  function. 

Equation  (8-1)  can  be  expressed  in  matrix  form: 
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af  0 

Zfs  0  0  Zfs  01  Z/s  0  2  •  Zfs  0k  •  Zfs  On 

Co 

O 

_ 1 

af  1 

Zfs  1  0  Zfs  1 1  Zfs  1  2  '  Zfs  1  k  '  Zfs  1  n 

Vs  l 

af  2 

Z fs  2  0  Zfs  21  Zfs  2  2  '  Zfs  2k  ’  Zfs  2n 

afv 

Zfs  p  0  Zfs  p  1  Zfs  p  2  Zfs  p  k  Zfs  p  n 

^5k 

afm_ 

Zfs  mO  Zfs  ml  Zfs  m2  '  Zfs  mk  '  Zfs  mn 

%  n 

(8-2) 


The  element,  ^  represents  submatrix  and  the  elements  a ,  ,  Sg  and  ^  subcol¬ 
umns  of  the  matrix  equation  (8—1).  The  submatrix,  Zy5  ^  represents  the  interaction  bet¬ 
ween  all  the  segments  /  and  s  in  the  antenna  for  a  fixed  basis  function,  k,  and  a  fixed 
boundary  condition  expansion  coefficient,  p.  The  subcolumn,  a y^,  represents  a  fixed 
boundary  condition  expansion  coefficient,  p,  for  all  the  segments,  /,  in  the  antenna.  The 
subcolumn,  6  q,  is  the  coefficient  of  the  triangular  current  basis  function  for  all  the  seg¬ 
ments,  5,  in  the  antenna.  The  subcolumn,  tjs  k,  represents  the  coefficient  for  a  fixed  charge 
basis  function,  (k  ^  0),  for  all  the  segments,  s,  in  the  antenna.  The  reason  for  this 
structure  will  be  discussed  later.  The  boundary  condition  is  imposed  by  setting  the  coeffi¬ 
cients,  OijQ  =  1/2  and  =  0  (m  ±  0),  for  each  segment  /. 

The  impedance  matrix,  Z jsmn,  depends  on  the  EMF(ty)  at  the  point  ty  on 
segment  /.  The  EMF(t y)  is  a  function  of  the  EMF  caused  by  all  the  current  and  charge 
distributions  on  all  segments,  s,  in  the  antenna: 


N-l  R  N-l  -R 

EMF(t ,)  =  E  £EMFsn(U)+  £  EEMF^JtA  (8-3) 

n=0  s  =  l  n  =  0  s  =  -1 


The  first  term  represents  the  EMF  of  one  arm  of  the  dipole.  The  second  term  is  the  sym¬ 
metric  image  term.  The  following  calculations  are  for  the  first  term  in  equation  (8-3). 
The  image  is  computed  with  the  same  method.  Note,  the  triangular  pulse  at  the  source  is 
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split  between  the  first  and  second  terms.  The  right  side  is  labeled  8^  and  the  left  side  is 


labeled 


th 


The  EMF  (t  f)  on  segment  /caused  by  the  n  basis  function  on  segment  s  is 
5  n  j 


EMFsn(t/)=4Si1(ti)-iSJl(‘/)-^ 


z(tf) 


(»-*) 


The  potential,  <j>s  (0),  is  the  potential  at  the  center  of  the  feed  point.  The  source  electric 
field  is  limited  to  a  gap  region;  if  the  gap  is  zero,  the  source  is  a  delta  function.  The  line 
integral  of  the  vector  potential  runs  from  the  center  of  the  feed  point  to  the  point  ty  on 
segment  /in  equation  (8-4). 

The  EMF  (t  «••)  in  equation  (8-4)  is  expanded  with  the  Legendre  polynomials  to 
give  the  impedance  matrix: 


M-l 

£  B  -EMF  (t,.). 
Mmi  s  nv  /  v 

i  =0 


(8-5) 


The  expansion  matrix,  /?  .,  is  computed  with  the  Mth  order  Guassian  integration  formula, 

equation  (7—7).  The  impedance  matrix  is  computed  by  substituting  equation  (8—4)  into 

equation  (8—5).  This  computation  of  the  impedance  matrix  naturally  divides  into  two 

parts:  the  scalar  and  vector  potential  terms.  Variables  attributed  to  the  scalar  and  vector 

potential  terms  in  the  impedance  matrix  are  identified  by  the  respective  superscripts,  s 

and  v.  The  scalar  and  vector  potentials  can  be  expressed  as  a  function  of  the  potential 

integral,  k  ,  developed  in  section  6.  A  different  notation  is  used  for  the  potential  inte- 
s  n 

gral,  Kg  n(ty-).  The  field  point  is  indicated  by  the  point  t y.  on  the  surface  of  segment  /. 
The  source  is  the  n^  basis  function  on  segment  s. 


69 


c 

The  scalar  potential  terms  in  Z  ^  mn  divide  into  two  cases:  the  triangular  basis 
function,  8  n,  and  the  Legendre  polynomial  basis  functions,  77  .  Each  triangular  current 

basis  function  contributes  two  square  pulses  to  the  scalar  potential  terms  in  submatrices 
Z/s  mO  *n’  e(llia^on  (8-2).  The  general  Legendre-polynomial  charge  basis  functions  are 

c 

used  to  compute  the  scalar  potential  terms  in  the  general  term,  Zj§  mn,  in  equation  (8—2). 

Each  current  basis  function  gives  two  square  pulses  for  the  charge  distribution.  The 
pulse  on  the  source  segment  s  has  a  height  —8  q/L.  The  other  pulse  is  on  the  previous 
segment,  s—1,  with  a  height  8  q/L.  Note,  for  the  source  triangular  current  basis  function, 
s  =  1,  the  q  term  is  included  in  the  image  term: 

M— 1  - 

E  B  •  -A- 

nmi  47ren 
i  =0  0 

_ g  g 

\Ks  Qpfi^s  +  0(t/i)-/c(5-l)  o(°^“TT.  5 

(8— 6  a) 


J/s  mO 


for  n  =  0  and 


M-l  , 

—  E  P  ■  -7 -  [k  (0)1 77 

'iru  47renL  s  nv  fv  s  nv  's  n 
i  =0  0  J 


for  n  >  0. 


The  notation  is  simplified  by  introducing  a  new  variable: 


(8— 6b) 


fs  mn 


M-l 
=  E 

i  =0 


@miKs 


(8-7) 


The  constant,  k  (0),  gives  only  one  expansion  coefficient,  'L 

o  XI 

pendent  of  the  boundary-condition  segment,  /: 


c 

fs  On 


This  constant  is  inde- 


fs  mn 


■WW0)' 
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The  ^mO  is  the  usual  Kr°necker  delta  Smm  =  1  and  tf  =  0  for  p  #  m.  In  this  new  nota¬ 
tion,  equation  (8-6)  reduces  to 


Z/smOPm(V  i?e0  _  ^f(s-l)  mO  ^  f (s-f)  m0_  E  ^/s  mO  ^/s  mO  U  x 

Pm(V)tf5  0  forn  =  0’  (8— 8  a) 

Z/s  mnPm^/^  ~  _  i^e0  ^ fs  mn“  ^  fs  mn  ^ s  nPm^/^  f°r  n  >  °' 

(8— 8b) 

The  factor,  Pm(t j),  is  included  in  the  equations  (8-8)  and  (8-9)  to  emphasize  that  mn 
is  the  pm(ty)  expansion  coefficient  of  the  EMF  on  segment  /.  The  notation  can  be  simpli¬ 
fied  by  replacing  the  quantity  in  the  square  bracket  with  a  single  variable,  T^  mn: 

Z/s  mnPm^/)  ~  ~~  i?e0T/s  mn^s  nPm(V^’  (8“9) 


The  vector  potential  terms  in  the  impedance  matrix  are  computed  by  numerically 
integrating  the  vector  potential.  The  numerical  evaluations  of  integrals  use  polynomial 
approximations  of  the  function.  Therefore,  it  would  be  natural  to  expand  the  integrand  in 
terms  of  Legendre  polynomials.  The  Legendre  polynomials  can  be  integrated  term  by  term 
to  give  Legendre  polynomials.  If  both  segments  / and  s  are  straight,  then  the  final  equation 


is  a  simple  function  of  ^ 


fs  mn' 


The  triangular  current  pulse  is  used  to  compute  the  first  column  of  the  impedance 


matrix.  The  higher  order  current  basis  functions  are  used  to  compute  the  remaining  col¬ 
umns  of  the  impedance  matrix.  The  current  pulse,  8  (segment  s  and  s-1 ),  contribute  the 
vector  potential  terms.  Note,  for  the  source  triangular  current  basis  function,  s  =  1,  the 


Ks—1  0  anc^  Ks—1  1 terms  are  included  in  the  image  term. 
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As°{tbi>  ~  iLhf 


«s- 


K  0  (‘a)  “  KS  iM**  +  [  Vi  o(‘«) +  Vi 1 M  Vi 

for  n  =  0.  (8— 10a) 


For  the  general  case  n  >  0,  the  vector  potential  on  segment  &  is  given  by  a  general  expres¬ 
sion, 


47rc2e1 


1  • 


2n+l 


'sn 


(8— 10b) 


The  line  integral  in  equation  (8-4)  is  evaluated  by  integrating  each  segment  individually: 


Z(tf) 


S-i 


ASR(z)-iz=Z 


i=l 


Ab  n(z)-d*  + 

s  e  gment  b 


‘*(t  f) 
AsjSz>iz 
*(-i f  ) 


(8-11) 


The  line  integrals  for  segments  1  to  segment  f-1  are  integrated  over  the  entire  segment. 
The  last  integral,  the  second  term,  is  integrated  over  part  of  the  segment. 

The  line  integral  element,  dz,  is  converted  to  u^dz  for  each  segment: 


’z(ty) 

ASnW-dz  =  lE, 


A o  (z)'U,dz  + 
s  nw  o 

s  e gment  & 


•*(* /)  , 

A$n(z)-ufAz. 

Z<-1  / ) 


(8-12) 


In  the  general  case,  the  integrand  Ag  n(z)  •  u ^  is  expanding  in  terms  of  Legendre  polyno¬ 
mials.  The  Legendre  polynomials  can  be  directly  integrated  by  using  equations  (4—6)  and 
(4-7). 
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If  the  unit  vectors,  u ^  ,  are  constant  on  each  segment,  the  A  vector  can  be  expanded 
with  Legendre  polynomials: 


*...«*»»=  rrs*1 nV‘6). 

47TC  e0  q  =  o 


(8-13) 


where 


1  1  1V1-JL 

fV^rhq„pq(‘S)=  SV.W'K))- 

47tc  e0  i  =o 


If  the  vector,  u  ,  implicit  in  the  vector  potential,  Ac  ( z ),  is  a  constant,  then  the  expansion 
s  s  n 


constant,  T^s  can  be  expressed  in  terms  of  s  ^ 


^ bs  qOPq^6  ^  bs  qO  ^ bs  ql^s  +  ^b(s-l)  qO  +  ^b(s-l)  ql 


K-JW 

(8-14) 


Tbs  qnPq^V  2n+1  bs  q(n-f-l)  ^bs  q(n-4)J  usPq^l)' 


(8—15) 


The  integral  of  the  vector  potential  is  simplified  by  substituting  equation  (8—13)  into  equa¬ 
tion  (8—12)  and  removing  the  constant  factors  from  the  integral.  The  dot  product, 
V  r^5  ,  can  be  removed  from  the  integral: 


*(*;,)  . 

Asn(z)-ubdz 

*H») 


*(-ij) 


z(‘6> 


(8-16) 
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For  a  segment  of  length,  L,  the  integration  of  Pq(xj)  can  be  expressed  as  the  product  of  a 

matrix,  S  ,  and  a  column  vector,  P  ft,): 

1  qm  ’  mv  bJ 


P  (2z/I)dz  =  4s  Pm(tj),  (8-17) 

J  m  =  0 


where 


_ 1_ 

2q+l 


m  m 


For  q  >  0 


=6™  +  81 


Om  ~u0 


•m 


For  q  =  0. 


The  integral  simplifies  further  if  t  ^  =  1: 

Pq(2z/L)dz  =  L(5q.  (8-18) 

The  integral  of  the  vector  potential,  A$  ,  from  the  feed  point  to  segment  /is  the  sum  of  / 
terms: 


■z(t  f) 

Asn(*)-iz  = 


1-1 


M 


,  V  Tbs  0.2Po(V)  +  S  V  Tfs  4"  S  VW 

1  q  =  0  in  =0 


(8—19) 
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The  terms  in  the  sum  from  b  =  1  to  /-  1  in  equation  (6—19)  are  constants;  they  contribute 
to  the  expansion  coefficient  of  Pg(t^ )  of  equation  (8—19).  The  vector  potential  terms  in 
the  impedance  matrix  is 


P  (t ,)  = 
fs  mn  m'1  f’ 


-1  w2L 


47rc2en  4  \  b=l 


f-U  1Vi  . 

S  u,-T\  n  2  6™  +  S  urT\  S  P  (tA 
b  bs  On  0  ^  /  fs  qn  qmj  mv  f ’ 


(8-20) 


The  terms  in  the  sum  from  b  =  1  to  /—  1  in  equations  (6—20)  and  (6—21)  only  contribute 
to  the  Pg(ty  )  expansion  term  in  the  impedance  matrix  (<5^  =  0  for  m  f  0). 

The  impedance  matrix  for  a  single  basis  function  on  one  segment  is 


where  n  <  N  and  m  <  N  for  N  basis  functions  and  ^boundary  conditions  per  segment. 
The  impedance  matrix  would  include  the  contribution  from  all  segments  on  both  arms  of 
the  dipole.  The  symmetry  in  the  dipole  can  be  used  to  eliminate  the  coefficients  of  the 
image  basis  function  as  a  variable,  5_s  =  8  and  r]_s  n  =  — (— l)n?/5  n  where  —  s  denotes  the 
image.  At  the  feed  point,  5=1,  the  current  pulse  is  split  between  the  image  and  the 
antenna.  At  free-wire  ends,  the  current  pulse  must  vanish. 


The  basis  functions  and  boundary  condition  discussed  in  this  report  allow  a  detailed 
analysis  of  the  errors  in  the  solution.  The  first  step  in  this  error  analysis  is  to  characterize 
the  scalar  potential  terms  in  the  impedance  matrix.  The  magnitude  of  each  term  is  deter¬ 


mined  by  the  expansion  coefficients,  ^ js  mn,  for  the  potential  integral,  k$  n(ty)-  Pirst  we 
will  examine  the  static  potential  integral,  k  (tA  contribution  to  the  magnitude  of  the 

S  u.  J 
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expansion  coefficient,  ^ js  mn-  The  plots  in  figures  6—2  and  6—3  show  that  the  potential  on 

the  surface  of  the  source  segment  has  the  same  shape  as  the  expansion  functions,  Pm(t  ); 

consequently,  the  self-term,  s$  nn,  for  the  charge  basis  function,  P  (t  ),  will  have  the 

largest  expansion  coefficient.  The  other  terms,  $  mn  m  =  n±2,  n±4,  etc.,  will  be  small 

for  the  self-term.  The  expansion  coefficients,  'P  =0,  are  for  the  even  and  odd 

!  ss  mn  1 

potentials  when  m  =  n±l,  n±3,  etc. 

q 

On  the  other  hand,  the  very  sharp  drop  in  the  static  potential  integral,  «“(t  ^  ),  cm 
the  adjacent  segment  contributes  to  all  of  the  expansion  coefficients,  ^s±-qs  mn,  for  all  m. 
The  plots  in  figures  6— 2b  and  6— 3b  show  that  the  scalar  potentials  and  their  expansion 
coefficients,  'P(5±(^5  mn,  are  much  smaller  for  d  >  2.  At  large  distances,  the  dipole  poten¬ 
tial,  -^(t y),  and  the  dipole  potential  created  by  each  current  basis  function  will  be  the 
dominate  terms  in  the  impedance  matrix.  The  frequency-dependent  part  of  the  potential 
integral,  «^(t j),  is  smooth  function  on  all  segments.  The  expansion  coefficients  of  the  fre¬ 
quency-dependent  part  of  the  potential  integral  term  contributes  terms,  ^  mn,  with 
m  =  0,  1,  and  2. 

The  computation  of  the  boundary  condition  ignores  expansion  terms,  $  r  mn,  with 
m  >  N.  For  the  self-term,  these  terms  will  be  small  when  the  proper  aspect  ratio  is  used  for 
each  segment  (table  7—1).  Any  errors  introduced  by  neglecting  these  terms  will  be  a  high- 
order  multipole  moment.  The  high-order  charge  basis  functions  in  the  impedance  matrix 
interact  primarily  with  nearby  segments. 

The  magnitude  of  the  current  terms  depend  on  the  constant,  w2L/4c2.  An  addi¬ 
tional  factor  of  L  from  the  vector  potential,  from  equation  (8—15)  (n  >  0),  can  be  included 
in  this  factor.  The  combined  constant  is  (wL/2c)2  for  n  >  1.  In  terms  of  the  wavelength, 
A,  the  constant  reduces  to  (rL/X)2.  The  contribution  of  the  high-order  current  basis  func¬ 
tions  to  the  impedance  matrix  can  be  made  smaller  by  decreasing  the  segment  size.  For  10 
segments  per  wavelength,  this  constant  is  1/10.  For  electrically  small  antennas,  the  cur- 
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rent  term  is  a  small  correction  to  the  static  part  of  the  impedance  matrix.  The  high-order 
basis  functions  model  the  near  fields  as  quasi-static  interaction. 

The  magnitude  of  the  current  term  is  determined  by  the  expansion  coefficients, 
fs  mn,  of  the  current  basis  functions.  The  expansion  coefficients  can  be  estimated  from 
the  current  potential  integral  plotted  in  figures  6-4  and  6—5.  The  potential  integral  com¬ 
puted  from  the  current  basis  functions  does  not  have  a  sharp  drop  at  the  end  of  the  source 
segment.  In  contrast  to  the  charge  basis  function,  the  potential  is  small  and  smooth  near 
the  end  of  the  segment.  On  the  surface  of  nearby  segments,  the  potential  is  a  smooth  func¬ 
tion;  it  can  be  described  with  a  few  low-order  Legendre  polynomials,  m  =  0,  1,  and  2. 
However,  on  the  surface  of  the  source  segment,  the  current  basis  function,  Jsn,  gives  two 
expansion  coefficients,  'P  ,  and  \P  ,  .  .  , .  Note:  the  correct  evaluation  of  the 

potentials  on  the  adjacent  segments  requires  both  potential  integrals,  ks  n_j.(t  )  and 
^sn+l^s)'  The  inteSrati°n  tw0  expansion  coefficients  give  three  terms: 

\P„„„  and  ^  ^  The  terms,  0  and  \Pon  „  .  0  „  ,  *>, 

ss  n— zn~ 2’  ss  nn’  ssn+2n+z  ’  ss  n—2  n—2  ss  n+zn+z’ 

interact  with  charge  basis  functions,  Pn _ 2(^5)  anc^  Pn+2^s^ 

In  the  calculation  of  the  impedance  matrix,  the  high-order  current  term  must  be 
approximated.  In  a  model  with  N  basis  functions  per  segment,  the  first  N  expansion  coeffi¬ 
cients,  ^ss  qq  to  ^5s  are  included  in  the  boundary  condition.  The  expansion 

coefficients,  ^P^  nN  to  ^5S  n+1  N+1’  introchiced  by  the  current  basis  functions, 
J$  N_2(ts)  and  Js  are  not  included  in  the  impedance  matrix.  This  approximation 

causes  an  PN  and  Pn-|-i  error  the  EMF  boundary  condition.  The  boundary  condition 
error  can  be  calculated  with  the  extra  expansion  coefficients,  'P ^  mn,  where  N  <  m  <  M. 

The  solution  is  computed  with  the  submatrix,  Z;  (m  <  N).  The  submatrix,  Z; 

c  js  mn  v  J  Js  mn 

(N  <  m  <  M),  is  used  to  compute  the  error  a ^  from  the  solution.  The  magnitude  of  the 
off-diagonal  elements  in  the  impedance  matrix  depends  on  the  constant,  (xL/2A)2. 

In  principle,  the  order  of  the  terms  in  the  impedance  matrix,  equation  (8—2),  does 
not  affect  the  accuracy  of  the  solution.  In  practice,  the  numerical  roundoff  errors  in  the 
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solution  process  can  create  large  errors  in  the  numerical  solution  (Press  et  al.,  1986).  The 
numerical  accuracy  of  the  solution  can  be  improved  using  partial  or  full  pivoting.  Pivoting 
moves  the  large  elements  to  the  diagonal  of  the  matrix.  Thus,  the  large  self-interaction 


terms,  Zss  ,  in  the  impedance  matrix  should  be  moved  to  the  diagonal.  However,  the 

self-interaction  terms,  Z  , ,  with  p  f  k  are  small:  these  terms  should  be  off  the  matrix 

diagonal.  Equation  (8—2)  has  the  large  submatrices,  Z  ,  ,  on  the  diagonal.  In  complex 

js  pp 

problems,  the  large  submatrices,  Z  f  ,  may  require  additional  pivoting.  For  a  dipole,  the 

Js  PP 

large  elements  are  already  on  the  diagonal  of  the  submatrix,  Z  ^  Gaussian  elimination 

is  used  on  the  small  upper  triangular  part  of  the  matrix.  Back  substituting  is  used  to  com¬ 
pute  the  final  solution. 


Figures  8—1  and  8—2  show  the  logarithm  of  the  magnitudes  of  the  impedance  matrix 

elements  for  a  dipole  at  resonances.  The  dipole  parameters  are  the  total  length,  2h  =  lm; 

the  wire  radius,  a  =  4.5401E— 5;  and  input  frequency,  146.0  MHz.  The  colors,  red,  purple, 

blue,  green,  yellow,  and  red  indicate  a  range  of  values  from  1  to  IE-6.  Matrix  elements 

smaller  than  IE— 6  are  plotted  as  black.  Figure  8—1  shows  a  56-segment  model  with  one 

basis  function  per  segment;  the  impedance  matrix  is  the  submatrix,  Z ^  qq.  Figure  8—2 

shows  an  8-segment  model  with  seven  basis  functions  per  segment.  In  figure  8—2,  45.7%  of 

the  matrix  elements  are  black.  In  Figure  8—2,  the  diagonal  elements,  Zgs  nn,  in  the  matrix 

are  self-interaction  terms  for  each  basis  function  on  each  segment,  the  4^  nn  expansion 

terms.  The  peaks  parallel  to  the  diagonal,  Zgs  n(n±2)  are  the  current  interaction  terms, 

4*  r  ;  and  the  self-interaction  terms,  4*  „  .  The  other  terms  in  the  peaks, 

ss  (n±2jn’  ’  ss  mn  ’ 

Z5(5±i)  nm,  are  the  nearest  neighbor  interactions,  4^^  nm.  The  peaks  in  the  impe¬ 
dance  matrix  model  the  near  field.  The  elements  between  the  peaks  model  the  far  field. 
The  far  field  is  modeled  by  the  Pq  and  basis  functions.  The  structure  of  the  impedance 
matrix,  in  figure  8—2,  should  improve  the  accuracy  of  the  solution  to  the  impedance 


matrix. 
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Figure  8-2.  A  logarithmic  plot  of  the  magnitude  of  the  elements  in  the  impedance  matrix  for  a 
8-segment  model  with  seven  basis  function  per  segment.  The  antenna  parameters  in  figure  2-4 
are  used  in  this  calculation. 


Insight  into  these  errors  can  be  gained  by  using  Singular  Value  Decomposition 
(Press  et  al.,  1986).  A  matrix  problem  can  be  solved  by  using  the  eigenvectors  and  eigen¬ 
values  of  the  matrix.  The  ith  eigenvalues,  A-,  and  eigenvectors,  x-,  obey  equation  (8—22): 

Aixi  =  Ax..  (8-22) 

In  this  discussion,  we  will  assume  that  A-  i  0.  The  solution,  S,  to  the  matrix  problem, 

Y  =  AS,  (8-23) 

can  be  found  by  expanding  column  vector,  Y,  in  a  sum  of  eigenvectors 

N— 1  . 

Y=  £  S:X.,  (8-24) 

i=0 

where 

s.  =  Y-x..  (8-25) 

The  solution  is 

N— 1  s. 

S=  S  -rix..  (8-26) 

i  =0  l  1 

The  equation  allows  the  effects  of  roundoff  errors  to  be  estimated.  The  calculation  of  the 
expansion  coefficient,  s-,  will  include  roundoff  error.  This  error  is  weighted  by  the  1/A-.  A 
small  eigenvector  will  cause  a  large  error  in  the  solution,  S.  The  numerical  stability  of  a 
matrix,  the  condition  number,  is  measured  by  the  ratio  of  the  largest  and  smallest  eigen¬ 
value.  The  condition  number  is  always  larger  than  one: 
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(8-27) 


max' 


mm 1 


A  high-condition  number  matrix  can  have  large  numerical  errors. 

The  condition  number  is  used  to  compare  the  numerical  stability  of  the  impedance 
matrices  plotted  in  figures  8—1  and  8—2.  The  condition  numbers  for  figures  8—1  and  8—2 
are  397.00  and  89.96,  respectively.  A  single-precision  impedance  matrix  is  used  to  compute 
the  condition  number.  The  high-order  basis  functions  and  boundary  condition  simplifies 
the  impedance  matrix  and  reduce  the  condition  number. 

In  general,  one  would  expect  the  condition  number  to  be  a  function  of  both  the 
number  of  segments  and  the  number  of  basis  functions  per  segment.  Table  8-1  shows  the 
condition  number  as  a  function  of  the  number  segments  (6,  12,  and  24)  and  of  the  number 
of  basis  functions  per  segment  (1  to  7).  The  existing  program  is  limited  to  97  unknowns. 


Table  8—1.  Condition  number  of  dipole  at  resonance. 


Number 
of  Basis 
Functions 

6  Segs 

1 

69.51 

2 

68.36 

3 

70.06 

4 

70.14 

5 

70.50 

6 

70.59 

7 

70.76 

12  Segs. 

24  Segs 

120.90 

208.23 

119.69 

206.3 

123.77 

215.66 

123.95 

124.79 

124.99 

125.37 
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Table  8—1  shows  that  the  condition  number  is  independent  of  the  number  of  basis  functions 
used  on  each  segment.  The  condition  number  depends  only  on  the  submatrix,  Z ^  qq.  The 
high-order  basis  functions  can  be  added  to  model  the  fine  details  in  the  charge  distribution 
without  increasing  the  condition  number. 

The  condition  number  for  a  dipole  at  antiresonance  is  much  smaller.  The  condition 
number  at  antiresonance  is  calculated  for  a  12-segment  model.  Table  8—2  shows  seven 
different  models. 


Table  8—2.  Condition  number  of  a  dipole  at  antiresonance. 

Number  of  12  Segs. 

Basis  Functions 


1 

2 

3 

4 

5 

6 
7 


6.3002 

6.273 

6.4906 

6.5077 

6.5566 

6.5728 

6.5965 


The  impedance  matrix  is  stable  at  antiresonance.  The  current  and  charge  distribution 
should  be  accurate  on  the  antenna.  However,  at  the  feed  point,  the  current  will  be  very 
small.  A  small  numerical  error  in  the  overall  solution  would  have  a  larger  impact  at  the 
feed  point.  The  numerical  results  in  section  2  show  that  this  is  not  a  problem  with  our 
algorithm. 


6.  D.  Fern,  personal  communication 
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The  Legendre-polynomial  basis  functions  are  closely  related  to  the  power  series  basis 
functions  used  by  Djordjevic  et  al.  (1990).  The  two  sets  of  basis  functions,  in  principle,  are 
equivalent;  however,  the  numerical  structure  of  the  impedance  matrices  has  important 
differences.  Each  power  series  basis  function  has  a  1/jr—  r'  |  term  in  its  potential  integral. 
The  net  charge  on  the  segment  determines  the  magnitude  of  the  1  / 1  r—r'  |  term.  For  a 
basis  function  tm  (0  <  t  <  1),  the  net  charge  is  l/(m+l).  The  details  modeled  by  the  high- 
order  power  series  basis  functions  are  lost  in  the  numerical  noise  at  large  distances.  The 
Legendre-polynomial  basis  function,  P  ,  replaces  the  1/|  r—r'  |  term  with  a  much  smaller 
1/|  r—r'  I11"*"*  term.  Each  Legendre-polynomial  basis  function  represents  a  distinguishable 
part  of  the  near  fields.  The  condition  number  is  independent  of  the  number  of  Legendre 
polynomial  basis  functions  used  on  each  segment. 

The  equal-width,  square-pulse  testing  boundary  condition  used  by  Djordjevic  et  al. 
(1990)  can  be  compared  to  the  EMF  boundary  condition  by  using  a  matrix  transformation 
of  the  impedance  matrix,  equation  (8—1).  The  calculation  of  this  matrix  transformation 
uses  two  steps.  The  impedance  matrix  is  used  to  estimate  the  EMF  on  each  segment.  The 
EMF(t)es£  on  segment  /is 


™^f)est 


M  N-l 
=  £  E 

m  =  0  n  =  0  S 


R 

E  Zr  P  (tr)x  , 
fs  mn  mv  fJAsn’ 


(8-28) 


where  M  is  the  number  of  Gaussian  quadrature  points  used  to  calculate  the  impedance 
matrix.  The  high-order  error  estimation  terms,  M  >  N,  improve  the  accuracy  of  the  EMF 
estimate. 

The  equal-width,  square-pulse  testing  boundary  condition  can  be  evaluated  from  the 

EMF(‘/W 
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(8-29) 


^).d,=  EMF(t/i+1)^-EMF(t/i)esi> 

z(k) 


where  t-  and  t-_^^  are  the  end  points  of  the  pulse  testing  function.  The  square-pulse  testing 

til 

boundary  condition,  BC^  for  the  i  pulse  on  the  segment  /is 


BC3?.  =  S  E  £  Z,  [P  (t*.  in)-P  (t,.)]x  • 
/ 1  _1  fs  mn  1  mv  fi+V  mv  /i'JAsn 

m  =  0  n  =  0  <S—  i 


M  N-l  R 


(8-30) 


The  new  impedance  matrix  for  the  point  match  boundary  condition  is  the  coefficient  of  x 


ZP,  •  =  E  Zr  fP  (tr-  M)-P  (t,.)l. 

fs  m  fs  mn  1  nr  fi+lJ 


(8—31) 


This  transformation  is  the  sum  of  the  rows  of  the  matrix  where  the  quantity  in  the 
square  bracket  is  the  weighting  coefficient  for  each  row.  The  effect  of  the  transformation  is 
to  add  the  large  diagonal  term,  Zgs  nn,  to  the  small  off-diagonal  terms  in  the  matrix.  The 
addition  of  the  large  term  causes  additional  roundoff  in  the  matrix  element.  To  solve  the 
matrix,  this  transformation  is  reversed;  rows  are  subtracted  to  eliminate  off-diagonal  ele¬ 
ments  in  the  matrix.  The  combined  effect  of  the  transformation  and  solution  of  the  matrix 
is  an  increase  in  roundoff  error  in  the  solution  of  the  impedance  matrix. 

The  equal-width,  square-pulse  testing  boundary  condition  will  not  work  for  high- 
order  Legendre-polynomial  basis  functions.  The  equal-width,  square-pulse  testing  functions 
are  wider  than  the  oscillation  of  k  n(t  )  near  the  end  of  the  segment.  The  EMF  boundary 
condition  can  resolved  the  fine  details  modeled  by  the  high-order  Legendre-polynomial  basis 
functions.  The  EMF  boundary  condition  isolates  the  interaction  of  the  basis  functions  on 
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different  segments.  The  stability  of  this  impedance  matrix  is  illustrated  in  the  condition 
number  calculations  and  numerical  calculations. 
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9  CONCLUSIONS 


The  error  in  the  method  of  moments  depends  on  the  models  used  to  approximate 
the  charge  distribution  and  the  boundary  condition.  The  modeling  errors  in  the  charge 
distribution  and  the  boundary  condition  interact  through  the  scalar  and  vector  potentials. 
The  basis  functions  and  the  boundary  condition  are  selected  to  minimize  the  interaction  of 
the  two  modeling  errors.  The  algorithm  used  to  impose  the  boundary  condition  also  esti¬ 
mates  the  error  in  the  boundary  condition. 

The  basis  functions  model  the  charge  distribution  on  each  segment  with  a  multipole 
expansion  of  the  charge  distribution.  The  first  basis  function  models  the  net  charge  on 
each  segment.  The  higher  order  basis  functions  are  used  to  redistribute  the  charge  to  elimi¬ 
nate  the  parallel  component  of  the  electric  field  on  the  segment.  The  higher  order  basis 
functions  model  the  detail  in  the  near  field  of  the  segment.  The  far  fields  of  the  higher 
order  basis  functions  are  small.  The  charge  distribution  on  one  segment  primarily  interacts 
with  the  nearby  segments. 

The  modeling  errors  in  the  boundary  condition  implies  an  error  in  the  charge  distri¬ 
bution.  The  charge  distribution  error  caused  by  the  boundary  condition  is  a  high-order 
multipole  moment.  The  boundary  condition  modeling  error  on  different  segments  are  iso¬ 
lated  from  each  other.  The  basis  functions  and  boundary  condition  modeling  errors  are 
local  to  each  segment.  This  creates  a  more  stable  impedance  matrix.  The  condition  num¬ 
ber  of  the  impedance  matrix  is  independent  of  the  number  of  basis  functions  per  segment. 
The  condition  number  depends  only  on  the  number  of  segments  in  the  model. 

The  computed  charge  distribution  agrees  with  the  theoretical  calculations  of  the 
charge  distribution.  The  high-order  basis  functions  give  the  best  model  of  the  charge  distri¬ 
bution  near  the  end  of  the  wire.  The  charge  distribution  converges  much  quicker  for  the 
high-order  basis  functions.  For  small  matrix  models,  the  high-order  basis  functions  are 
more  accurate. 
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For  large  matrix  models,  the  convergence  of  the  admittance  depends  on  the  matrix 
size  used  to  model  the  dipole.  The  current  at  the  feed  point  of  the  dipole  depends  only  on 
the  height  of  the  triangular  current  basis  function.  The  height  is  computed  very  accurately 
by  the  boundary  condition  algorithm.  Even  the  square  and  linear  charge  basis  function 
models  are  accurate  at  resonance  and  antiresonance.  This  shows  that  the  error  in  the 
boundary  condition  is  a  high-order  multipole  moment. 

The  numerical  examples  discussed  in  this  report  have  high  aspect-ratio  segments. 
For  lower  aspect-ratio  segments,  the  shape  of  the  potentials  near  the  end  of  the  segment 
impacts  the  evaluation  of  the  boundary  condition.  Additional  work  needs  to  be  performed 
to  understand  the  limitation  introduced  by  lower  aspect-ratio  segments.  The  current  work 
could  be  extended  to  model  general  symmetric  dipoles  without  loops. 

The  basis  functions  and  boundary  condition  simplify  the  structure  of  the  impedance 
matrix.  The  far  field  interaction  is  modeled  by  the  triangular  current  basis  functions.  The 
near  field  is  modeled  by  the  high-order  basis  functions.  This  approach  combines  the  flexi¬ 
bility  of  the  full-wave  approach  with  the  stability  of  the  quasi-static  approach.  The  algo¬ 
rithm  is  accurate  for  both  the  near  and  far  fields  of  the  antenna.  The  error  estimation 
should  be  a  useful  tool  in  evaluating  the  errors  in  the  electric-field  boundary  condition. 


88 


10  REFERENCES 


Butler,  C.  M.  1975.  "Evaluation  of  Potential  Integral  at  Singularity  of  Exact  Kernel  in 
Thin-Wire  Calculations,"  IEEE  Transactions  on  Antennas  and  Propagation,  AP— 23(Mar), 
pp.  293-295. 

Canning,  F.  1992.  "Sparse  Approximation  for  Solving  Integral  Equations  with  Oscillatory 
Kernels,"  Siam  Journal  Scientific  and  Statistical  Computing,  vol.  13,  no.  l(Jan),  pp.  71—87. 

Djordjevic,  A.  R.,  M.  B.  Bazdar,  G.  M.  Vitosevic,  T.  K.  Sarkar,  and  R.F.  Harrington. 
1990.  Analysis  of  Wire  Antennas  and  Scatterers:  Software  and  User’s  Manual.  Artech 
House,  Boston,  MA. 

Dudley,  D.  G.  1995.  "Error  Minimization  and  Convergence  in  Numerical  Methods," 
Electromagnetic,  vol.  5,  pp.  89—97. 

Feynman,  R.  P.,  R.  B.  Leighton,  and  M.  Sands.  1964.  JThe  Feynman  Lectures  on  Physics, 
Volume  II.  pp.  18—10,  Addison-Wesley  Publishing  Co.,  Reading,  MA. 

Gradshteyn,  I.  S.  and  I.  M.  Ryzhik.  1980.  Table  of  Integrals,  Series,  and  Products.  Cor¬ 
rected  and  Enlarged  Edition.  Academic  Press,  Inc.  Orlando,  FL. 

Harrington,  R.  F.  1968.  Field  Computation  by  Method  of  Moments.  Macmillan  Company, 
New  York. 


89 


Hwu,  S.  U.  and  D.  R.  Wilton.  1988.  "Electromagnetic  Scattering  and  Radiation  by  Arbi¬ 
trary  Configurations  of  Conducting  Bodies  and  Wires,"  NOSC  TD  1325  (Aug),  Naval  Com- 

* 

mand,  Control  and  Ocean  Surveillance  Center  RDT&E  Division,  San  Diego,  CA. 

Jackson,  J.  D.  1975.  Classical  Electrodynamics,  John  Wiley  &  Sons,  Inc.,  New  York,  NY. 

King,  R.  W.  P.  1956.  The  Theory  of  Linear  Antennas;  with  Charts  and  Tables  for  Practical 
Applications,  Harvard  University  Press,  Cambridge,  MA. 

Logan,  J.  C.  and  J.  W.  Rockway.  1986.  "The  New  MININEC  (Version  3):  A  Mini-Numeri¬ 
cal  Electromagnetic  Code,"  NOSC  TD  938  (Sep),  Naval  Command,  Control  and  Ocean 

* 

Surveillance  Center  RDT&E  Division,  San  Diego,  CA. 

Michalski,  K.  A.  and  C.  M.  Butler.  1983.  "An  Efficient  Technique  for  Solving  the  Wire 
Integral  Equation  with  Non-uniform  Sampling,"  IEEE  Southeastcon  Conference 
Proceedings  (April). 

Popovic,  B.  D.,  M.B.  Dragovic,  and  A.  R.  Djordjevic.  1982.  Analysis  and  Synthesis  of  Wire 
Antennas.  Research  Studies  Press,  New  York,  NY. 

Press,  W.  H.,  B.  P.  Flannery,  S.  A.  Teukolsky,  and  W.  T.  Vetterling.  1986.  Numerical  Rec¬ 
ipes  the  Art  of  Scientific  Computing.  Cambridge  University  Press,  Cambridge,  United 
Kingdom. 


*  Formerly  Naval  Ocean  Systems  Center  (NOSC). 


90 


Richmond,  J.  H.  1965  "Digital  Computer  Solutions  of  the  Rigorous  Equations  for  Scatter¬ 
ing  Problems,"  Proceedings  of  the  IEEE  (pp.  796—805).  August. 


Watt,  A.  D.  1967.  VLF  Radio  Engineering,  pp.  403,  Pergamon  Press,  New  York,  NY. 


APPENDIX  A:  NUMERICAL  CALCULATION  OF  POTENTIAL  INTEGRAL 


In  equation  (6—2),  the  z'  integration  can  be  simplified  by  integrating  the  potential 
integral  by  parts  n  times  in  z'  and  using  equation  (4—8).  As  a  simple  example,  consider 

«J(p,z): 


p2  7T  rL  j  2 
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«l(Az)  =  2i  WM- - dz'd(P'i 

Jo  J-L/2  lr_r,l 


(A— 1) 


where  |  r—r'  |  =  [(z— z')2+p2+a2— ! 2apcos^/]'*'^.  Integrating  by  parts: 
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The  first  term  is  zero  at  ±  L/2: 
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The  derivative  simplifies  to 
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where  cos0  = - .  The  angle  6  is  measured  from  the  z-axis. 
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At  large  distance,  |r|  >>  L/2,  0  and  j  r-r'  \  are  almost  constant, 
approximated  as 
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for  | r |  >>  L/2. 


This  integral  is  the  dipole  potential  at  large  distances. 
In  the  general  case, 
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The  derivative  of  the  static  potential  integral  can  be  reduced  to 
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The  integral  can  be 
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The  final  integral  is 


1 

•27T 

r1 

'(l-z'2)nl 

Pn(cos(0)) 

2? r 

0  * 

1 

to 

i _ 

|r-r'jn+1 

dz'd^'. 


(A— 12) 


The  far-field  form  of  the  potential  integral  is  explicit  in  equation  (A— 12).  At  large  dis¬ 
tances,  the  integrand  is  a  slowly  varying  function;  the  numerical  integration  is  very  accur¬ 
ate.  Equation  (A— 12)  is  also  used  to  calculate  the  potential  integral  close  to  the  wire. 

In  the  near  field,  the  equation  (6—2)  must  be  directly  evaluated.  In  this  case,  the 
term  ln((l+r)/(l— r))  is  not  accurate  for  r  near  1.  Several  alternate  forms  exist  for 
ln((l+r)/(l— r))  for  r  >  .5: 
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and 
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where 
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Equation  (A— 15)  can  be  simplified  by  separating  the  numerator  and  denominator: 
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The  first  term  is  almost  constant;  the  numerical  evaluation  converges  very  quickly.  The 
second  term  can  be  analytically  evaluated: 
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The  integral  is  the  potential  integral  of  an  infinite  cylinder  of  radius  a.  This  integral  can 
be  evaluated  with 
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where  c  =  p2+ a2  and  d  =  2a p.  The  square  root  term  simplifies  to  v  c2-d2  =  <J(p2— a2)2. 
The  square  root  is  p2- a2  for  p  >  a  and  a 2-p2  for  p  <  a: 
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for  p  <  a. 
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This  is  the  expected  result  for  an  infinite  cylinder.  The  term,  Ln[p2+a2— 2apx.],  in  equa¬ 
tion  (A— 16)  can  be  replaced  by  Ln[p2]  for  p  >  a  or  Ln[a2]  for  p  <  a.  Equation  (A-15)  redu¬ 
ces  to 
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Equations  (A— 20a)  and  (A— 20b)  are  valid  for  |z|  <  L/2;  however,  when  |z|  =  L/2 
and  p  =  a,  the  numerator  vanishes.  In  this  case,  equations  (A-13)  and  (A-14)  are  used  to 

compute  the  potential.  For  z  =  L/2,  the  term,  R2,  reduces  to  7p2+a2-2apxi.  For 

z  =  —L/2,  the  term,  R^,  reduces  to  Jp2+a2—2apxy. 
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The  procedure  used  to  simplify  equations  (A— 20)  gives 
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For  the  case,  p  >  a, 
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Equations  (A— 13),  (A— 14),  (A— 20),  (A— 21),  and  (A— 22)  give  the  best  convergence  and 
accuracy  in  the  region  close  to  the  wire. 

In  the  near  field,  the  potential  integral  can  be  computed  in  closed  form.  The  exact 
calculation  involves  large  cancellations  at  moderate  distances  from  the  wire.  Only  minor 
simplifications  are  possible. 
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